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Abstract 


The  extent  of  buried  munitions  and  explosives  of  concern  (MEC),  or  unexploded  ordnance 
(UXO),  creates  a  serious  environmental  hazard  in  the  U.S.,  especially  where  military  bases  and 
ranges  are  being  converted  to  civilian  use  (millions  of  acres).  The  cleanup  of  these  sites  presents 
a  formidable  challenge  to  DoD,  most  prominently  in  terms  of  accurate  and  reliable  detection  and 
classification.  Discrimination  between  MEC  and  relatively  safe  background  clutter  (scrap  metal) 
is  a  monumental  problem  that  depends  on  the  instrument  technology,  as  well  as  on  the 
processing  methodology  that  inverts  the  detection  data  to  infer  MEC.  Reducing  the  number  of 
false  alarms  in  the  detection  process  can  save  billions  of  dollars  in  the  cleanup  effort.  That 
process  and  the  consequent  savings  rely  in  large  part  on  precise  three-dimensional  positioning 
(geolocation)  of  the  detection  sensors  both  to  aid  in  the  inversion  of  the  data  in  post-survey  data 
processing  and  in  the  mapping  and  recovery  of  locations  with  positive  identification.  Although 
the  detection  instrument  technology  has  advanced  significantly  in  the  last  decade,  it  is  often  still 
the  geolocation  technology  that  defines  or  limits  the  accuracy  of  MEC  detection. 

The  objective  of  this  project  was  to  develop  and  test  novel  geolocation  algorithms  applied  to 
scenarios  typical  of  MEC  detection  and  recovery,  where  the  precision  goals  are  1  cm  and  10  cm 
for  three-dimensional  positioning  of  magnetic  and  electro-magnetic  detection  devices.  The 
principal  geolocation  system  was  assumed  to  be  based  on  inertial  measurement  units  (IMUs) 
integrated  with  differential  GPS.  The  GPS  receiver  presumably  is  of  geodetic  quality 
commensurate  with  the  precision  goals  and  serves  as  the  reference  system.  The  IMUs  are 
tactical-grade,  three-dimensional  accelerometers  and  gyroscopes  that  supplement  GPS 
positioning  by  a)  improving  the  temporal  (hence  spatial)  resolution;  and  b)  bridging  temporary 
GPS  outages  or  degradations  caused  by  obstructions  of  the  satellite  lines-of-sight  or  caused  by 
other  electronic  interferences,  as  well  as  multi-path  effects.  A  major  component  is  also  the 
calibration  of  systematic  sensor  errors  using  external  measurements  and  information.  The 
resulting  techniques  and  software  were  tested  on  ground  vehicle  and  man-portable  systems. 

Alternatives  to  the  standard  extended  Kalman  filter  were  developed  and  tested.  Using 
simulated  data,  as  well  as  data  obtained  from  actual  instrumented  systems  in  the  laboratory  and 
in  the  field,  we  tested  various  optimal  non-linear  filters  and  smoothers  in  order  to  demonstrate 
the  interpolation  capabilities  of  medium-grade  IMUs.  The  unscented  Kalman  filter  performed 
significantly  better  than  the  standard  version,  particularly  over  highly  dynamic  curved  segments 
of  the  simulated  and  actual  trajectories,  yielding  up  to  50%  improvement  in  the  position 
accuracy.  Similar  improvement  was  obtained  for  the  unscented  particle  filter,  and  its  adaptive 
variant,  over  the  unscented  Kalman  filter  when  the  statistical  distribution  of  the  IMU  noise  was 
non-symmetric  (i.e.,  essentially  non-Gaussian).  While  the  few-centimeter  geolocation  accuracy 
goal  for  highly  dynamic  UXO  characterization  applications  remains  a  challenge  if  tactical  grade 
IMUs  are  integrated  with  a  significantly  degraded  ranging  system,  using  filters  appropriate  to  the 
inherent  nonlinear  dynamics  and  potential  non-Gaussian  nature  of  the  sensor  noise  tend  to  reduce 
overall  errors  compared  to  the  traditional  filter. 
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1.  Introduction 


1.1  Background 

In  North  America,  especially  in  the  US,  the  problem  of  unexploded  ordnance  (UXO)  extends 
beyond  military  training  sites,  fields  of  battle,  and  military  bases.  According  to  the  Defense 
Science  Board  Task  Force  on  UXO,  15  million  acres  or  more  are  potentially  contaminated  at 
over  1500  different  sites  (Defense  Science  Board,  2003).  The  typical  method  to  find  UXO 
employs  time  domain  electromagnetic  induction  (EMI)  devices  that  transmit  pulsed  magnetic 
fields  and  then  sense  the  fields  produced  by  induced  currents  in  buried  highly  conducting  objects 
such  as  UXO.  However,  the  false  alarm  rate  is  often  very  high  and  the  detection  and  remediation 
of  UXO  in  contaminated  sites  using  such  detect- and-dig  methods  is  extremely  inefficient  since  it 
is  difficult  to  distinguish  the  buried  UXO  from  the  noise  of  geologic  magnetic  sources  or 
anthropic  clutter  items  such  as  exploded  ordnance  fragments  and  agriculture  or  industrial 
artifacts  (Bell,  2001).  The  cost  can  be  as  high  as  $1.4  million/acre  (Collins  et  al.,  2001;  Zhang  et 
al.,  2003).  Thus,  new  aiding  technologies  that  can  improve  discrimination  performance  for  UXO 
detection  are  urgently  needed. 

Digital  Geophysical  Mapping  (DGM)  integrates  multiple  sensors  to  detect  and  characterize 
potential  UXO.  However,  without  highly  precise  position  data  DGM  is  ineffective  because  the 
position  uncertainty  can  degrade  the  inversion  of  EMI  data  and  the  subsequent  classification 
(Tarokh  et  al.,  2004;  U.S.  Army  Corps  of  Engineers  2006).  The  requirements  of  position 
accuracy  for  UXO  detection  fall  into  three  levels;  initial  screening  that  needs  tens  of  centimeters 
accuracy  as  standard  deviation,  area  mapping  that  needs  less  than  5cm,  and  characterization  and 
discrimination  that  need  less  than  2cm  accuracy. 

The  Global  Positioning  System  (GPS)  has  been  used  as  a  precise  positioning  system 
specifically  for  geolocation  of  UXO  detectors.  Laser-based  ranging,  radio,  and  acoustic  ranging 
have  also  been  considered  as  alternative  positioning  methods.  All  ranging  systems  perform 
geolocation  by  trilateration,  which  requires  unobstructed  lines-of-sight  and  uninterrupted 
transmission  of  signals  that  connect  the  detector  and  at  least  three  reference  points.  From  this 
standpoint,  GPS  is  representative  of  all  ranging  systems.  However,  the  temporal  resolution  of 
GPS  kinematic  positioning  (1  Hz,  typically)  is  less  than  for  some  other  systems,  such  as  laser 
ranging,  which  diminishes  the  efficiency  of  UXO  detection  system  that  pass  over  target  sites  at 
high  speed  (Bell,  2005). 

The  problem  of  bridging  gaps  in  the  ranging  solution  due  to  obstructed  or  otherwise 
interrupted  signals  is  solved  by  integrating  GPS  (or  any  other  ranging  system)  with  an  inertial 
measurement  unit  (IMU),  which  is  an  autonomous  relative  positioning  device.  Indeed,  a  recently 
conducted  comprehensive  program  to  assess  the  geolocation  technology  for  UXO  detection  and 
remediation  (U.S.  Army  Corps  of  Engineers,  2006)  concluded  that  all  geolocation  systems  could 
benefit  from  augmentation  with  an  IMU.  Since  IMU’s  generally  also  provide  very  high  temporal 
resolution  (50  Hz  -  250  Hz),  and  since  GPS  is  the  most  efficient  and  cost-effective  ranging  tool, 
the  integrated  IMU/GPS  ranging  system  is  the  target  of  this  study  for  uninterrupted,  high- 
resolution,  high-accuracy  geolocation  in  support  of  UXO  detection  and  characterization. 


1.2  Technical  Objectives 

The  three  main  data  processing  steps  of  an  integrated  IMU/ranging  system  are  data 
preprocessing,  data  integration,  and  post -processing  with  smoothing.  This  project  focuses  on  all 
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three  steps,  where  the  first  includes  de-noising  algorithms  and  self  calibration  of  IMU  systematic 
errors  in  the  field,  and  the  remaining  steps  are  investigated  from  the  viewpoint  of  optimal 
filtering  and  smoothing.  The  latter  include  the  (extended)  Kalman  filter,  and  various  nonlinear 
filters,  such  as  the  unscented  Kalman  filter  and  particle  filter,  as  well  as  investigations  of  the 
adaptive  Kalman  filter  and  the  introduction  of  neural  networks  to  improve  the  calibration  of 
errors  and  the  overall  geolocation  solution. 

Developing  optimal  algorithms  to  integrate  the  data  from  ranging  systems  and  inertial 
measurement  units  has  occupied  the  navigation  community  for  many  years;  and,  the  (Extended) 
Kalman  Filter  (EKF)  has  been  the  workhorse  for  IMU/GPS  integration  for  several  decades. 
However,  new  filters  recently  introduced  and  tested  have  the  potential  for  significant 
improvements  (Aggarwal  et  al.,  2006)  over  conventional  methods  that  must  contend  with  linear 
approximations  of  intrinsically  nonlinear  dynamics.  These  new  filters,  particularly  the  unscented 
Kalman  filter,  are  the  subject  of  the  present  investigation,  applied  to  the  precise  geolocation  of 
ground-based  MEC  detector  platforms.  Since  the  trajectories  of  MEC  detectors  can  be  highly 
dynamic  with  rapid  changes  in  orientation  and  acceleration,  filters  that  obviate  any  assumption  of 
linearity  should  perform  better.  In  addition,  particle  filters  have  recently  been  implemented  that 
avoid  the  usual  assumption  of  Gaussianity  in  the  sensor  errors  (Haykin,  2001;  Gustafsson  et  al., 
2001;  Arulampalam  et  ah,  2001;  Simon,  2006).  We  consider  also  adaptations  of  these  filters  that 
allow  for  unknown  variations  in  the  sensor  and  observation  noise  statistics. 

Other  studies  that  have  addressed  the  potential  improvement  of  the  unscented  Kalman  filter 
over  the  traditional  extended  Kalman  filter  in  navigation  applications  include  Van  der  Merwe 
and  Wan  (2004)  who  predicted  up  to  30%  improvement  on  an  unmanned  airborne  vehicle 
(UAV).  St.  Pierre  and  Gingras  (2004)  and  Shin  (2005)  showed  only  slightly  better  or  mixed 
results  based  on  simulations  for  a  land  vehicle  IMU  analysis.  Yi  and  Grejner-Brzezinska  (2006) 
also  obtained  improvement  when  GPS  is  blocked  (free-inertial  navigation),  but  there  was  no 
improvement  during  periods  of  GPS  coverage. 

However,  any  of  these  filters  cannot  overcome  the  natural  accumulation  of  process  errors  as 
the  temporal  gap  in  ranging  solutions  increases,  especially  in  dynamic  environments.  One 
method  to  reduce  the  divergence  between  IMU-determined  and  true  trajectories  is  to  apply  a 
post-processing  smoothing  between  GPS  updates,  such  as  the  Rauch-Tung-Striebel  (RTS) 
backward  smoother,  which  is  based  on  the  extended  Kalman  filter  (EKF).  We  have  tailored  this 
type  of  smoothing  to  the  nonlinear  filters  (Lee  et  al.  2008),  and  for  this  project  the  resulting 
adaptive  EKF-RTS  smoothing  (AEKS)  and  adaptive  UKF-RTS  smoothing  (AUKS)  methods 
were  tested  and  analyzed  using  data  from  an  inertial  geolocation  system.  In  addition,  we  tested  a 
simple  end-matching  algorithm  in  place  of  the  smoothing,  which  portends  further  studies  in 
modeling  position  errors  under  in  real  time  that  approaches  the  post-processing  precision  with 
known  constraints  on  the  dynamics  of  the  vehicle. 

We  also  developed  and  tested  a  neural  network  as  an  aid  to  both  the  adaptive  extended 
Kalman  filter/smoother  (NN-AEKS)  and  the  adaptive  unscented  Kalman  filter/smoother  (NN- 
AUKS)  to  bridge  outages  in  the  ranging  solution  for  UXO  detectors  in  quiescent  and  moderately 
dynamic  environments.  Finally,  we  considered  the  use  of  dual  IMUs  in  the  integrated  system 
with  the  aim  of  self-calibrating  non-common  IMU  errors  in  the  post-processing  mode  using  a 
wave  correlation  filter  (WCF). 
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2.  Technology  Description 


2.1  The  Kalman  Filter 

The  Kalman  filter  is  an  optimal  (minimum  mean  square  error),  linear,  recursive  estimation 
method  and  is  amply  treated  in  the  literature,  as  is  its  extension  to  deal  with  non-linear  models 
(Haykin,  2001;  Kalman,  1960;  Brown  and  Hwang  1992;  Grewal  and  Andrews  1993).  Its 
application  to  the  IMU/GPS  integration  problem  is  also  well  documented  (Greenspan,  1996; 
Farrell  and  Barth,  1998;  Rogers,  2000;  Jekeli,  2000;  Titterton  and  Weston,  2004).  The  extended 
Kalman  filter  (EKF)  is  founded  on  two  main  assumptions:  linearity  in  (or  linearization  of)  the 
system  dynamics  and  observation  models  and  Gaussian  noise  excitation  (Gordon  et  al.,  1993; 
Gelb,  1974;  Grewal  and  Andrews  2001).  These  two  factors  bring  several  limitations  and 
implementation  difficulties  to  the  integration  problem,  such  as  the  derivation  of  the  Jacobian 
matrices  for  complicated  non-linear  systems  and  potentially  unrealistic  statistical  models  of  the 
errors  (Julier  et  al.,  1995).  Also,  only  relatively  small  error  states  are  allowed  in  the  EKF, 
otherwise  the  first-order  approximations  can  lead  to  instabilities  in  the  form  of  biased  solutions 
and  inconsistency  in  the  covariance  updates. 

The  basic  linear  model  of  a  system  of  error  states,  xk ,  at  epoch,  tk ,  is 


^  (4 ’4-1  )  Xk- 1  +  Gk  Wk  , 

(1) 

Hkxk  +  vk , 

(2) 

where  <4>  is  the  state  transition  matrix;  Gk  is  an  appropriate  scaling  factor  for  the  driving  zero- 
mean  Gaussian  white  noise,  wk ;  and  yk  is  an  observation  linearly  related  through  Hk  to  the 
states  and  also  includes  zero-mean  Gaussian  white  noise,  vk .  The  error  states  typically  include 

position  errors,  velocity  errors,  and  orientation  errors  associated  with  the  navigation  solution 
derived  from  IMU  data,  and  are  augmented  by  various  systematic  errors  assigned  to  the  IMU 
sensors.  The  observation  (errors)  are  associated  with  external  navigation  aids,  such  as  GPS, 
which  provide  information  that  allows  for  the  estimation  of  the  error  states  of  the  system. 

The  Kalman  filter  has  two  steps:  first,  a  prediction  of  the  states  at  epoch,  tk ,  given  their 

values  at  epoch,  tk_x ;  and,  second,  an  update  at  tk  based  on  the  observation  at  that  epoch.  In  the 
prediction,  the  state  estimate  and  its  covariance  propagate  according  to: 

=<^  (4’4-i)**-i’  (3) 

Pk  =<J)  {tk,tk_i) Pk-{&  (4 ’4-i )  +  GkQkGk  ,  (4) 

where  Qk  is  the  covariance  of  wk ,  Pk_x  is  the  covariance  matrix  of  the  states  at  the  epoch,  tk_ x . 
Updating  the  states  on  account  of  the  observations  yields 

<5> 
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Pk=(l-KkHk)Pk(I-KkHk)r  +  KkRkK [  (6) 

where  Rk  is  the  covariance  of  the  observation  noise  and  the  Kalman  gain  is  given  by 

Kl=Pl-Hl(HkPkHl+Rl)~'  (7) 

For  later  use,  we  also  define  the  innovation  and  its  estimate,  as  in  equation  (5), 

V  =  A  -  Hkxk  (8) 

W  =  yk  -  Hkx~k  •  (9) 


2.2.  Sampling-Based  Filters 

The  sampling-based  filtering  methods  differ  from  the  standard  Kalman  filter  in  that  they 
propagate  sample  states  rather  than  their  statistics  through  the  system  model.  As  such,  no 
linearization  is  required,  nor  is  it  required  to  formulate  and  compute  the  attendant  derivatives  of 
the  system  model.  Depending  on  the  statistical  assumptions  for  the  states,  the  propagated 
samples  are  analyzed  to  infer  the  propagated  statistics.  We  consider  two  sampling-based  filters, 
the  unscented  Kalman  filter  (UKF,  also  called  sigma-point  Kalman  filter)  and  the  particle  filter 
(PF).  The  UKF  avoids  the  linearization  associated  with  the  EKF,  but  still  assumes  Gaussian 
state  variables,  and  the  estimation  of  the  propagated  statistics  is  achieved  on  the  basis  of  the 
unscented  transformation  (Wan  and  Van  Der  Merwe,  2001).  The  PF  obviates  both  the 
linearization  and  the  Gaussianity  assumption  by  propagating  a  sufficient  number  of  samples 
based  on  the  Monte  Carlo  approach  to  statistics  estimation. 


2.2.1  Unscented  Kalman  Filter  (UKF) 

The  development  of  the  unscented  Kalman  Filter  was  initiated  by  Julier  and  Uhlmann  (1996)  and 
Julier  et  al.  (1995;  2000).  The  UKF  is  a  recursive  application  of  the  unscented  transformation 
(UT)  of  state  variables  through  the  nonlinear  state  dynamic  model.  The  UT  propagates  a  suitably 
chosen  set  of  sample  points  (called  sigma  points)  in  the  state  space  through  the  (nonlinear) 
system  dynamics  such  that  they  accurately  capture  the  transformed  mean  and  covariance  of  the 
states.  We  use  the  more  general  scaled  version  of  the  UT  given  by  Wan  and  Van  der  Merwe 
(2001).  For  the  random  variable  x  (dimension,  nx)  with  mean,  x,  and  covariance,  Px,  the 

2 nx  + 1  sigma  points  are  generated  as  follows 
%o=* 

=*+a(V(n*+K)jP0.,  (io> 

Xi=x-a^(nx+K)Px^,  i  =  nx  +  l,...,2nx 
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where  a  and  k  are  scaling  parameters  and  +K )  Px  j  is  the  /'Lh  row  or  column  of  the  matrix 

square  root  of  (nx  +K )  Px  .  Given  a  nonlinear  function,  g  (x) ,  it  can  be  shown  that  the  following 
weighted  combinations  of  y,=g(%,)  estimate  the  first  two  statistical  moments  (mean  and 
covariance)  of  g  at  least  up  to  second  order  in  the  non-linearities: 

2 nx  2  nx 

y=Lw<">yi’  p,=L,w^^-y)^-y)T ■  <u> 

i-  0  i= 0 

where  the  corresponding  weights  are  given  by 


W0H  =  1- 


a2{nx+K) 


,  W0{c)=l 


ty.(m)  -  = 


2a2  ( nx  +  k  ) 


«2k+K) 

,  i  =  l,...,2nx 


+  ^  1  —  cx. 2  + 13  j 


(12) 


The  sum  of  the  weights  for  the  mean  is  unity,  while  for  the  covariance  they  sum  to  (l  -a2  +  |3 ) . 

The  scaling  parameter,  a  (10“4  <a  <  1),  controls  the  spread  of  the  sigma  points  around  x  and 
serves  to  maintain  the  positive  semi-definiteness  of  the  covariance  (Wan  and  Van  Der  Merwe, 
2001).  For  Gaussian  states,  x,  the  estimation  of  the  mean  and  covariance  of  g  is  accurate  up  to 
third  order.  The  scaling  parameter,  ,  in  addition  to  a ,  is  used  to  increase  the  accuracy  of 
higher-order  moments  (for  Gaussian  variables,  J3  =  2  is  optimal).  The  scaling  parameter,  k  , 
was  set  to  K  =nx  in  our  simulations,  and  we  varied  a  . 

The  linear  propagation  and  observation  of  states  defined  by  equations  (1)  and  (2)  are  replaced 
by  the  more  general  non-linear  forms;  however,  we  still  assume  additive  noise: 

xk=f{xk.l)  +  Gkwk,  (13) 


yk=h{xk)  +  vk-  (I4) 

The  Kalman  filter  using  the  unscented  transformation,  i.e.,  the  UKF,  then  proceeds  with  the 
usual  two-step,  prediction  and  filter  formalism,  but  the  mean  and  covariance  are  determined  from 
the  sigma  points. 


2.2.2  Unscented  Particle  Filter  (UPF) 

The  particle  filter  was  introduced  already  in  the  1950s,  but  was  forgotten  due  to  the  lack  of 
computing  power  (Godsill  et  al.,  2000).  However,  with  modem  computational  capabilities,  the 
particle  filter  has  also  been  applied  to  INS/ranging-system  integration  (Angermann  et  al.,  2006) 
and  to  the  problem  of  transfer  alignment  (Hao  et  al.,  2006).  The  performance  of  the  PF  generally 
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improves  as  the  number  of  sample  states  (particles)  used  to  capture  their  probability  density 
function  (pdf)  increases,  which,  however,  also  increases  the  computational  burden.  There  are 
three  types  of  particle  filter,  the  generic  (or  traditional)  particle  filter,  the  unscented  particle  filter 
(UPF,  also  called  sigma-point  particle  filter),  and  the  extended  particle  filter.  We  considered 
only  the  UPF. 

Particle  filters,  also  known  as  bootstrap  filtering,  the  condensation  algorithm,  Monte  Carlo 
filtering,  and  survival  of  the  fittest  methods,  implement  Bayesian  estimation,  which  under  the 
special  case  of  Gaussian  noise  and  linear  models  yields  the  familiar  Kalman  filter  (Gordon  et  al., 
1993).  The  object  is  to  determine  the  probability  density  of  the  state  at  time,  tk ,  conditioned  on 
the  measurements  up  to  that  time,  according  to  Bayes’  Rule  (Simon,  2006,  p.464) 


p{xkhk,yk-vyk-2’---)  = 


p{yk  \yk-i,yk-2’---) 


(15) 


where  the  conditional  density,  p ( yk  lx,.),  is  presumed  known  (for  example,  but  not  necessarily 

Gaussian),  and  the  other  conditional  densities  on  the  right  side,  in  principle,  can  be  determined 
from  the  previous  steps  in  the  recursive  algorithm. 

In  practice,  the  particle  filter  starts  with  the  randomly  generated  particles,  x0i  ,i  =  1,...,  N , 
using  their  a  priori  mean  and  covariance  and  propagates  them  through  equation  (13): 

xk,i  =f(xk-u  )  +  Gkwk,i’  (16) 


where  wk  i  is  the  noise  realized  according  to  its  known  pdf.  From  the  measurement,  yk ,  and  its 

known  pdf,  the  likelihood,  q,  =  p(yk  / xk_u),  can  be  computed  for  each  particle.  These  relative 

likelihood  values,  normalized  to  sum  to  unity,  are  then  used  to  obtain  a  re-sampling  of  the 
particles,  whose  pdf  (now  conditioned  on  the  observations)  tends  to  the  desired  Bayesian  density 
(equation  (15)).  From  this  the  mean  and  covariance  can  be  computed  in  the  usual  ways.  There 
are  several  re-sampling  techniques  and  particle  modifications  available  to  improve  the 
performance  of  the  filter.  We  used  the  sequential  importance  sampling  (SIS))  method  (Haykin, 
2001;  Wan  and  Van  Der  Merwe  2001). 

Van  der  Merwe  et  al.  (2000)  combined  the  PF  with  other  filters  such  as  UKF  in  order  to 
refine  the  initially  generated  particles,  and  called  it  the  unscented  particle  filter  (UPF).  Each 
randomly  generated  initial  particle  is  propagated  and  updated  by  using  the  UKF  (where  the 
sigma  points  are  computed  using  the  particle).  Then,  the  likelihood  values  are  determined  using 
these  a  posteriori  particles  and  the  measurement  pdf,  as  before,  followed  by  the  re-sampling 
procedure  (Haykin,  2001;  Simon,  2006).  This  approach  can  yield  a  better  approximation  for  the 
conditional  pdf  of  the  states  and  has  yielded  improved  accuracy  for  various  positioning-related 
applications  (Arulampalam  et  al.,  2001;  Grewal  and  Andrews,  2001;  Van  der  Merwe  et  al., 
2000). 
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2.2.3  Rao-Blackwellized  Unscented  Particle  Filter  (RBUPF) 

To  overcome  the  computational  inefficiency  of  the  PF,  the  states  can  be  divided  into  two  subsets 
(linear  and  non-linear)  and  this  type  of  filter  is  called  the  Rao-Blackwellized  particle  filter 
(RBPF)  (Nordlund  and  Gustafsson,  2001).  The  RBPF  applies  the  Kalman  Filter  to  the  linear  and 
Gaussian  states  of  the  system  model.  The  remaining  state  variables,  which  suffer  from  severe 
non-linear  and/or  non-Gaussian  structure,  are  solved  using  the  unscented  particle  filter. 

Partitioning  the  state  vector  xk  into  non-linear,  x\ ,  and  linear,  x2k ,  parts  the  system  and 

measurement  models  become  (compare  equations  (1)  and  (2)  for  the  linear  case  and  equations 
(13)  and  (14)  for  the  general  non-linear  case): 

xL  =  fl(4)+Fk4+wL 

xL=f2(xl)  +  Fk2x2k+wl  t  (17) 

yk=h(4)+Hk-xl+vk 


The  system  noise  vector  is  divided  into  two  noise  vectors,  w\  and  w2k  that  are  uncorrelated,  and 

have  known  covariance  matrices,  Qlk  =  e[w\  {w\ )  j  and  Ql  =  E |vt^  [w2k )  j ,  respectively.  To 

compute  the  a  posteriori  pdf  of  the  filter,  we  estimate  p(xk  I  yk )  =  p[x\,x2k  I  yk )  by  factoring  it 
into  two  parts  according  to  Doucet  et  al.  (2001)  and  Nordlund  (2002), 

p  (4 .  ■ xl 1  yk )  =  p  (4 1 4  ’  yk )  p  (4 1  yk ) .  ( 1 8) 

and  by  first  computing  each  pdf  on  the  right  side  according  to  the  appropriate  filter  methods. 


2.3  Adaptive  Filtering 

The  drawback  of  the  navigation  filters  discussed  so  far  is  that  they  are  adversely  affected  by 
potentially  inaccurate  system  and  measurement  noise  statistics.  Particularly  for  highly  dynamic 
trajectories  (such  as  may  be  encountered  by  UXO  detection  equipment)  the  error  statistics  may 
vary  in  time,  thus  invalidating  the  initially  defined  models.  To  improve  the  filters  under  such 
circumstances  and  mitigate  such  effects,  adaptive  methods  are  often  employed. 

Among  the  various  parameters  that  require  a  priori  specification  in  the  UKF,  including  the 
initial  states  and  their  covariance,  the  unscented  transformation  parameters,  and  the  covariances 
of  the  process  noise  ( Q )  and  measurement  noise  (R),  the  latter  influence  most  significantly  the 
filter’s  performance  and  stability.  The  initial  states  and  covariance  have  negligible  influence  as 
the  filter  processes  more  and  more  data.  The  UT  parameters  only  affect  the  higher  order  terms  of 
the  nonlinear  model  and  have  little  impact  on  the  accuracy  of  the  estimated  position  (Lee  and 
Jekeli,  2009).  Therefore,  usually  an  innovation  based  covariance  matching  algorithm  is 
employed  in  order  to  tune  the  Q  and  R  matrices  of  the  Kalman  filter. 

This  technique  is  based  on  the  supposition  that  the  actual  filter  residuals  should  be  consistent 
with  their  theoretical  covariance  (Jwo  and  Huang,  2004).  In  this  paper,  since  a  loosely-coupled 
INS/GPS  system  was  employed,  we  considered  only  the  system  (IMU)  noise.  According  to 
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Salychev  (1999),  the  a  priori  information  of  the  system  (as  represented  by  the  covariance 
matrices,  Q)  can  be  adjusted  according  to  the  accuracy  of  estimation. 

If  covariance  matrix,  Rk ,  is  unknown  or  inaccurate,  it  can  be  estimated  from  the  covariance 

matrix  of  the  innovation  sequence  (  X)k  ,  equation  (9)).  From  equation  (8),  with 
Pk=E[xk{x-k)T),  let 


Ck=E(vyk)  =  HP~HT+Rk, 


(19) 


Then,  let  an  estimate  of  this  covariance  be 


Q=-  V  ^>7’ 

m  .  “  , 

j=k-m+l 


(20) 


where  m  is  the  estimation  window  size,  which  satisfies  the  recursion: 


A  l^^r 

Q-  —  ]  Ck_ i  +  — vkvk  . 
k  k 

Substituted  into  equation  (19),  it  yields  an  estimate  of  the  measurement  covariance  matrix: 
Rk=Ck-HP~HT. 


(21) 


(22) 


In  case  the  covariance  of  the  system  noise,  Qk  ,  is  unknown,  it  can  be  shown  similarly  that 

GQtGT=KtE(vtvTt)KTt,  (23) 

from  which  we  approximate 

GkQkGTk  =  KkCkKTk  .  (24) 

Originally  developed  for  the  KF  (and  EKF)  this  technique  has  also  been  employed  for  the 
UKF  (Song  and  Han  2008).  We  have  also  applied  this  technique  to  the  particle  filter  and  call  it 
the  Adaptive  Unscented  Particle  Filter  (AUPF).  It  is  implemented  and  validated  by  the  tests  of 
the  simulations  described  in  Section  3. 


2.3.1  Neural  Network  Aided  Adaptive  Filtering 

The  adaptive  estimation  of  the  covariance  using  equation  (20)  depends  on  the  proper  choice  of 
the  window  size,  which  for  optimal  results  may  be  variable  and  difficult  to  determine  (Jwo  and 
Huang,  2004).  Because  of  this  potential  shortcoming,  we  considered  modifying  these  adaptive 
filters  using  methods  applied  in  artificial  intelligence  (AI).  AI  provides  a  successful  and 
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effective  solution  to  certain  engineering  and  science  problems  which  cannot  be  solved  by  using 
conventional  methods  (Cawsey,  1998).  Various  approaches  such  as  neural  networks,  fuzzy  logic, 
evolutionary  computing,  probabilistic  computing,  expert  program,  and  genetic  programming  can 
provide  an  “artificial  intelligence”  (defined  as  the  ability  to  learn,  understand  and  adapt)  to 
complex  and  uncertain  systems  (Honavar  and  Uhr,  1994).  Among  these  techniques,  we  selected 
the  neural  network  to  aid  adaptive  filters  because  it  can  learn  input-output  relationships  without  a 
priori  knowledge  of  the  dynamic  models  and  noise  statistics  of  the  measurements.  Employing 
the  neural  network  for  similar  applications  has  already  been  demonstrated  with  success  by  Jwo 
and  Huang  (2004),  Wang  et  al.  (2006,  2007),  Korniyenko  et  al.  (2005),  and  Chiang  and  El- 
Sheimy  (2004).  Specifically,  Wang  et  al.  (2006)  proposed  a  neural  network  and  Kalman  filter 
hybrid  approach  to  reduce  KF  drift  during  GPS  outages.  Also,  Wang  et  al.  (2007)  used  the 
neural-network-aided  adaptive  KF  to  reduce  vehicle  dynamic  variations  and  to  improve  the 
navigation  solution.  And,  Zhan  and  Wan  (2006)  derived  a  multi-layer,  neural-network-based 
unscented  Kalman  filter  for  nonlinear  estimation,  and  their  simulation  results  showed  overall 
improvement  of  the  filter  performance. 

A  neural  network  consists  of  neurons  that  process  an  input  to  generate  an  output  (Figure  1). 
The  neuron  includes  a  synaptic  weight  (or  weight),  an  adder,  and  a  transfer  (activation)  function 
(Haykin,  1999).  An  input  signal,  xi  (i  =  1,2, ...,«),  is  multiplied  by  a  synaptic  weight,  vv/( ,  and 

connected  to  neuron,  j .  An  adder  sums  up  the  weighted  input  signals  (simply,  v;  =  wjjxl  ). 
The  neuron  model  also  has  a  bias,  bj  (also  called  the  external  threshold)  that  is  used  to  increase 
or  decrease  the  input  to  the  transfer  function.  Therefore,  the  simple  model  of  a  neuron  is 


yj  =  f 


W.X: 
i  Jl  1 


i= 1 


(25) 


where  /  is  a  transfer  function  and  y  ;  represent  the  output. 

A  transfer  function  modifies  and  limits  the  amplitude  range  of  the  output  signal,  typically 
normalizing  it  in  the  range  of  [0,1]  or  [-1,1].  Common  transfer  functions,  selected  according  to 
the  application,  include  the  hard  limiter  (a  binary  or  bipolar  output),  a  linear  (or  piecewise  linear) 
function,  and  the  sigmoid  (s-shaped)  nonlinear  function  (Haykin,  1999).  We  use  the  sigmoid 

function,  f  (v)  =  (l  +  e~av}  ,  with  slope  parameter,  a  ,  because  the  relationship  between  input 

and  output  in  our  case  is  essentially  non-linear  (see  also  Ham  and  Kostanic,  2001). 
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Neuron 


The  neurons  are  used  to  form  layers  that  build  a  neural  network.  The  architecture  of  the 
network  can  be  classified  as  single-layer  feed-forward  (one  hidden  layer  and  an  output  layer; 
SFN),  as  multi-layer  feed-forward  (multiple  hidden  layers  and  an  output  layer;  MFN),  and  as 
recurrent  (at  least  one  feedback  loop  based  on  a  SFN  or  a  MFN;  RN)  (Haykin,  1999).  Demuth 
and  Beale  (2004)  showed  that  the  two-layer  feed-forward  network  (first  layer  is  sigmoid  and 
second  is  linear)  can  be  trained  to  approximate  any  arbitrary  nonlinear  function.  Golden  (1996) 
also  showed  that  an  MFN  can  be  designed  to  provide  the  best  approximation  accuracy  to  the 
unknown  model.  We  applied  the  multi-layered  feed-forward  neural  network  (three  layers)  to  aid 
the  various  adaptive  filtering  techniques  for  the  INS/GPS  system. 

The  usual  way  to  decide  on  the  appropriate  number  of  hidden  neurons  is  empirical  (however, 
see  also  Bishop,  1995;  and  Haykin,  1999).  That  is,  many  candidate  networks  having  different 
numbers  of  hidden  neurons  should  be  tested  to  determine  the  one  with  the  best  performance. 
Our  laboratory  tests  indicate  that  optimal  geolocation  results  are  obtained  with  16,  24,  and  15 
neurons,  respectively,  in  the  three  layers,  where  sigmoid  transfer  functions  are  used  in  the  first 
and  second  layers  and  the  third  layer  is  linear. 

The  principal  idea  for  this  project  is  the  hybrid  method  of  traditional  adaptive  covariance 
estimation  aided  by  a  neural  network  trained  on  given  platform  dynamics.  The  Kalman  filter 
estimates  the  navigation  errors  in  position,  velocity  and  attitude  using  external  control.  At  the 
same  time,  the  neural  network  is  trained  to  map  a  relationship  between  the  platform  dynamics 
(the  input)  and  the  Kalman  filter  estimations  (the  desired  output).  The  input  comprises  changes 
in  velocity  ( Avw,Av£,AvD )  in  a  local  north-east-down  coordinate  frame,  and  Euler  angles, 
(|),0,t|/  ,  for  the  platform  attitude,  and  their  changes,  A(|),A0,A\|/  .  The  Euler  angles  are 
determined  from  the  gyro  data  and  the  changes  in  these  and  in  the  velocity  are  calculated  from 
the  last  measurement  update  to  the  current  measurement  update: 


A vN  —  vN  k+l  vNk ,  Ave  —  vEk+l  vEk ,  Avn  —  vDk+l  vD  k 


A<l>  =<!>** H>*.  A0=0*+1-0*.  Av=v 


k+ 1 


-Yu 


(26) 


where  k  is  the  measurement  update  index. 
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Wang  et  al.  (2006,  2007)  suggested  that  rapid  changes  in  the  heading  angle  can  disturb  the 
training  of  the  neural  network.  However,  in  our  laboratory  tests,  training  with  the  heading  angle 
produced  better  results  because  this  angle  clearly  identifies  the  dynamic  maneuvering  of  the 
platform.  For  the  neural  network’s  desired  outputs  ( d ),  we  selected  the  innovations  of  the 
Kalman  filter,  given  by  equation  (9).  For  the  unscented  Kalman  filter  the  desired  output  of  the 
neural  network  is  yk—yk  ( yk  is  from  equation  (11)). 

If  measurements  (external  control  points)  are  available,  the  neural  network  is  trained  at  the 
control  update  rate  using  all  available  input  and  desired  output  values  until  it  reaches  a  certain 
pre-defined  error  (threshold).  The  weights  and  biases  of  the  network  were  adjusted  iteratively 
using  the  Levenberg-Marquardt  algorithm  to  minimize  the  differences  between  the  computed 
output,  y,  and  the  desired  output,  d  (Chiang  and  El-Sheimy,  2004;  Wang  et  al.,  2006).  When 
measurements  are  not  available  (during  a  GPS  outage),  the  computed  output  of  the  neural 
network  is  used  to  determine  the  process  noise  covariance  ( GQGT )  according  to  equations  (24) 
and  (20).  In  our  Laboratory  tests,  since  the  neural  network  does  not  have  enough  training  data 
from  the  first  few  control  points,  the  process  noise  ( Q )  is  estimated  by  the  traditional  adaptive 
filtering  method  (equation  (24))  until  the  neural  network  is  successfully  trained. 


2.4  Smoothing 

The  filters  described  above  are  recursive  algorithms  based  on  the  conditional  expectation  of  the 
state  given  all  observations  and  states  up  to  the  current  time  step  k.  In  contrast,  smoothing 
estimates  the  states  by  using  also  available  observations  in  the  future. 


2.4.1  Kalman  Smoother  (KS) 

From  given  observations  over  the  interval  0  <  k  <  N  for  fixed  N,  if  the  forward  and  backward 
estimate  ( x[  and  xhk )  and  their  error  covariances  ( Pk  and  Pk  )  are  available,  the  smoothed 

estimate,  xsk ,  and  its  covariance,  Pk  ,  are  obtained.  Since  it  is  assumed  that  the  process  noise 
wk  and  measurement  noise  vk  are  independent,  we  may  formulate  the  smoothed  estimate  and  its 
a  posteriori  error  covariance  matrix  as 


p;= 


=  rt([PtfT4+[itj' 


[*r+[*r 


X,. 


-1 


(27) 


This  smoother  (three-part  smoother)  has  three  components;  a  forward  filter,  a  backward  filter, 
and  a  separate  smoother  which  combines  results  embodied  in  the  forward  and  backward  filters 
(Haykin,  2001).  However,  The  Rauch-Tung-Striebel  (RTS)  smoother  differs  from  this  smoother 
in  that  the  measurements  are  processed  by  the  forward  filter  and  then  a  separate  backward 
smoothing  pass  is  used  to  obtain  the  smoothing  solution  (Sarkka,  2008).  Also,  the  Rauch-Tung- 
Striebel  smoother  is  more  efficient  than  the  three-part  smoother  in  that  a  single  entity  can 
perform  smoothing  by  incorporating  the  backward  filter  and  separate  smoother  (Rauch  et  al., 
1965). 
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2.4.2  Particle  Filter  Smoothing 

In  a  literature  review  on  particle  smoothing  or  Monte  Carlo  Smoothing,  it  is  (interestingly) 
noticed  that  although  particle  filter  theory  and  applications  are  frequently  treated,  there  are  only 
few  published  papers  on  particle  filter  smoothing  because  the  smoothing  algorithms  such  as  the 
two-filter  smoother  (TFS),  the  forward-backward  smoother  (FBS),  and  a  maximum  a  posteriori 
(MAP)  smoother  typically  incur  high  computational  costs,  namely  of  the  order  of  o(iV2), 

compared  to  the  PF  which  has  costs  of  the  order  of  0(n)  (Klass  et  al.,  2006). 

The  basic  idea  of  particle  smoothing  is  that  the  particle  filter  can  provide  a  smoothed  result 
automatically  if  the  whole  history  of  each  particle  of  states  is  stored  (Kitagawa,  1996).  That  is, 
from  the  filtered  particles  of  the  UPF,  each  particle  is  smoothed  by  using  the  unscented  RTS 

(Sarkka,  2006).  Each  smoothed  particle  can  be  defined  as  xf0 ,  i  =  l,...,N  (N is  the  number  of 

particles)  and  then  the  mean  value  of  the  smoothed  estimate  can  be  found  from  a  simple  average 
(Kitagawa,  1996;  Sarkka,  2008): 


~^knean  yy 


E 


^(0 


(28) 


2.4.3  Rao-Blackwellized  Particle  Smoother  (RBPS) 

For  the  smoothing  version  of  the  RBPF  we  can  employ  the  same  method  as  for  the  UPS,  which 
employed  the  unscented  RTS  algorithm  to  smooth  each  particle.  The  RTS  smoother  was  applied 
to  each  of  the  mean  and  covariance  histories  of  the  particles,  xj:'(, ,  P^ ,  i  =  1, . . . ,  N  to  produce  the 

smoothed  mean  and  covariance  x,.N ,  PtsN  . 


2.5  Wave  Correlation  Filter 

The  Wave  Correlation  filter  (WCF)  ideally  extracts  common  spectral  components  from  two 
signals,  while  rejecting  disparate  ones,  thus  removing  errors  if  the  two  signals  are  supposed  to 
refer  to  the  same  source.  Of  course,  common  errors  are  not  removed.  With  this  filter  we 
anticipate  obtaining  a  more  accurate  positioning  solution  from  the  dual  IMU  measurements, 
where  random  errors  left  in  the  two  final  solutions  should  be  uncorrelated  and  amenable  to 
elimination  or  reduction  by  the  WCF  (a  similar  approach  was  used,  e.g.,  by  Serpas  (2003)  and  by 
Li  (2009)  to  improve  solutions  derived  from  independent  data  sets).  Let  the  two  EKF  solutions 
along  the  trajectory  be  denoted  x,  (k)  and  x2  (k) ,  respectively.  The  corresponding  Fourier 

transforms  are  G,(l)  and  G2(l),  where  l  is  the  wave  number.  The  correlation  coefficient,  rn  per 
wave  number  is  obtained  by 


Re(G1(/))Re(G2(/))  +  Im(G1(/))lm(G2(/)) 
0,(1)  G2(i) 


(29) 
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Each  solution  is  filtered  according  to  a  pre-defined  threshold: 


r  /,\_JGu(0  r,^  threshold 

1,2  [  0  if  ri<  threshold 

The  final  position  estimates  are  obtained  by  transforming  the  average  of  the  retained  frequency 
components  back  into  the  time  domain: 

p(k)  =  ^F~l{Gl(l)  +  G2(l)},  (31) 

where  is  F  the  Fourier  transform. 


2.6  End-Matching 

The  remaining  trend  and  bias  errors  of  both  the  EKF  and  UKF  solutions  after  WCF  processing 
could  be  corrected  in  an  ad  hoc  manner  by  the  end- matching  method  with  respect  to  the  GPS 
positions.  In  terms  of  latitude,  longitude  and  height  coordinates,  the  end-matching  algorithm  is  a 
simple  linear  fit  to  given  data  at  the  ends  of  a  free- inertial  trajectory: 


(r 

K 

d) 

T  w 

K 

"i  (t-t0) 
l  (t-t0  ) 

X 

K 

_K_ 

_i  (t-t(>)_ 

m 

(32) 


where  <j \,Xr,hr  are  GPS  points,  §w,\w,hw  are  the  position  after  filtering  with  wave  correlation 
filter.  b0  is  a  bias  and  m  is  the  trend  as  a  function  of  epoch.  The  values  of  bias  and  trend 
(matrix  X)  can  be  found  by 

X  =(AT  A)~l  ■  ATY .  (33) 
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3.  Filter  Simulations 


The  initial  investigations  were  based  on  simulations  to  ascertain  the  general  performance  of  the 
different  filters  and  smoothers  under  typical  MEC  survey  geometries  and  dynamics.  Without 
loss  in  generality,  we  selected  the  loosely  coupled  INS/GPS  integration  scheme  based  on  the 
decentralized  filter  architecture  (Jekeli,  2000).  A  tightly  integrated  system  typically  would  be 
used  in  GPS  embedded  inertial  navigation  systems,  where  both  the  IMUs  and  the  GPS  aid  each 
other  to  obtain  an  optimal  blended  solution.  Looking  toward  our  eventual  applications  where  an 
isolated  IMU  would  be  added  to  an  existing  geolocation  system,  the  loosely  coupled  integration 
seemed  more  appropriate. 


3.1  Simulation  Setup 

Aside  from  the  different  filtering  methods,  a  number  of  system  and  environmental  factors  were 
considered,  including  the  IMU  sensor  quality,  the  dynamics  of  the  detector  platform  (typical  for 
either  hand-held  or  cart-mounted  deployment),  and  the  ranging  solution  quality  (concerning  both 
precision  and  environmental  degradation,  as  well  as  longer  outages).  The  IMU  sensor  may  be 
categorized  as  commercial  grade  (e.g.,  Crossbow  IMU400C),  tactical  or  medium  grade  (e.g., 
Honeywell  HG1700),  and  navigation  grade  (e.g.,  Litton/Northrop-Grumman  LN100).  Table  1 
offers  representative  specifications,  as  published  by  the  corresponding  vendors,  for  the  most 
important  errors  associated  with  the  accelerometers  and  gyroscopes  for  each  grade  of  inertial 
unit. 


Table  1.  IMU  levels  of  accuracy. 


IMU 

Accelerometer 

Gyroscope 

Bias 

Scale 

Random  walk 

Bias 

Scale 

Random  walk 

IMU400C 

8.5  mg 

104  ppm 

5  mg/  a [hr 

3600  °/hr 

104  ppm 

0.85  °/y[hr 

HG1700 

1.0  mg 

300  ppm 

0.25  mg  /  4hr 

1  °/hr 

150  ppm 

0.125  °/y[hr 

LN-100 

20  pg 

40  ppm 

5.0  ng/yfhr 

0.01  °/h 

1  ppm 

0.001  °/4hr 

The  two  levels  of  ranging  precision  are  associated  with  either  radar  or  laser  ranging, 
exemplified  by  real-time  kinematic  (RTK),  geodetic  quality,  differential  GPS  (e.g.,  Trimble  4700 
receiver)  and  the  geodetic  total  station  (e.g.,  Leica  TPS  1100  total  station).  Each  level  of 
precision  is  affected  by  the  distance,  respectively,  to  the  fixed  GPS  base  station  and  to  the  total 
station  (Clynch,  2001;  Grejner-Brzezinska,  2001;  Kim,  2004).  Table  2  lists  typical  precision 
levels  considered  in  our  analysis.  Further  degradation  in  the  ranging  solution  will  occur 
primarily  because  of  signal  occlusions  (although  other  factors,  such  as  tropospheric  and 
ionospheric  delays,  and  the  geometric  configuration  of  the  transmitters,  also  affect  the  signal  and 
solution  quality).  Signal  outages  may  occur  because  of  intervening  manmade  structures  or 
natural  objects  (usually  a  foliage  canopy  in  the  case  of  GPS;  also  rugged  terrain  in  case  of  the 
total  station).  We  assume  relatively  short-lived  outages  on  the  order  of  several  to  tens  of 
seconds.  It  should  be  emphasized  that  while  our  analyses  address  primarily  the  inertial  sensor 
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capability  to  aid  the  ranging  solution  in  an  operational  setting,  for  this  investigation  the  ranging 
solution  is  also  viewed  as  the  means  to  estimate  inertial  sensor  errors. 


Table  2.  The  assumed  precision  of  ranging  solutions. 


Performance  (s  r) 

RTK-DGPS 

Horizontal:  10mm  +  lppm 
Vertical:  20mm  +  lppm 

Total  Station 

Horizontal  &  Vertical: 

2mm  +  2ppm 

Four  different  filtering  algorithms  were  implemented  and  compared  in  our  simulations:  the 
extended  Kalman  filter  (EKF),  the  unscented  Kalman  filter  (UKF),  the  unscented  particle  filter 
(UPF),  and  the  adaptive  unscented  particle  filter  (AUPF).  In  order  to  assess  these  filters  in 
typical  dynamic  environments  for  MEC  detection  and  characterization,  the  accuracy  of  the 
integrated  positioning  system  was  determined  separately  for  curved  and  straight  sections  of  the 
trajectory.  We  also  considered  tuning  the  parameters  of  the  unscented  transformation,  and 
evaluated  the  benefits  of  using  the  particle  filters  when  the  input  noise  is  non-Gaussian. 

The  state  vector  for  the  IMU/ranging  integrated  solution  comprises  21  states:  three  position 
(latitude,  longitude,  height)  errors,  ( 5(| ),8X,8/j  );  three  velocity  errors,  (  8<| >,8X,8/i  );  three 
orientation  errors  in  a  local  north-east-down  frame,  (\|/ ^  ,\|/ ^  ,\|/ ^ );  and  a  bias  and  scale  factor 

error  for  each  of  the  three  accelerometers  and  three  gyros.  For  the  UKF  and  UPF,  instead  of  the 
three  orientation  angle  errors,  the  four  corresponding  quaternions  were  employed  because  the  UT 
operates  best  in  a  Euclidean  space  (i.e.,  this  avoids  the  transcendental  functions)  (Shin,  2005; 
Kraft,  2003). 

For  a  given  trajectory,  error-free  IMU  data  were  generated  using  the  Matlab  INSToolkit®. 
The  input  and  output  specifications  of  INSToolkit  are  described  in  (Jekeli  and  Lee,  2007).  The 
trajectory  was  defined  for  a  platform  following  a  planar  path  with  0.5  m/s  velocity  and 
meandering  sweeps  across  a  given  area.  The  generated  accelerometer  data  are  velocity 
increments,  Av ,  and  the  gyro  data  are  angle  increments,  A0  .  The  temporal  resolution  of  this 
control  trajectory  was  0.02  s  . 

The  actual  error-corrupted  IMU  data  were  then  simulated  using  errors  specified  in  Table  1 
according  to  the  models 

A0  =  (l  +  k(j )  A0  +  [d  +  Wq  )  At 

(34) 

Av  =  (l  +  kA )  Av  +  (b  +  wA )  At 

where  d  ,  kG,  and  wG  are  the  bias,  scale  factor  error,  and  random  noise  for  the  gyro;  and  b ,  kA, 
and  wA ,  are  corresponding  quantities  for  the  accelerometer.  The  time  interval  between  the 
generated  data  is  At  =  tk+l  -  tk  (0.02  s  in  our  simulations).  Corresponding  navigation  equations 

for  the  IMUs  were  integrated  using  these  corrupted  data  to  generate  the  trajectory  coordinates  as 
indicated  by  an  inertial  navigation  system. 
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In  addition,  the  ranging  solution  for  the  trajectory  was  determined  by  simply  adding  random 
noise  to  the  true  coordinates.  That  is,  we  assumed  the  ranging  system  to  be  free  of  all  systematic 
errors  and  the  final  solution  corrupted  only  by  white  noise: 

+v-ar 

K  =  Xr  +v  cr  (35) 

hr  =hr  +v  ar 

where  V  is  a  random  variable  following  the  zero-mean,  unit-variance  Gaussian  distribution,  and 
<5r  is  given  in  Table  2.  The  ranging  solution  coordinates  are  used  in  the  filter  to  update  the  IMU 
position  errors  and  to  estimate  the  systematic  errors  in  the  inertial  sensors.  The  Figure  2  shows 
the  simulated  trajectory,  designed  for  typical  MEC  mapping  and  detection,  with  18  180°  turns 
and  19  straight  segments.  The  turns  represent  highly  non-linear  dynamics  as  modeled  by  the 
navigation  equations  and  generally  pose  a  challenge  to  linearized  filters  such  as  the  EKF.  The 
ranging  updates  in  the  filter  can  be  introduced  at  longer  intervals  to  simulate  prolonged 
interruptions.  Figure  3  shows  the  general  simulation  and  analysis  process. 
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Figure.  2.  The  generated  control  path. 
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Figure  3.  Flow  chart  of  loosely  coupled,  decentralized  ranging/LMU  integration  and  simulation 
analysis. 


3.2  Results  and  Analysis 

The  free-inertial  navigation  solutions  obtained  just  prior  to  incorporation  of  the  external  updates 
using  the  EKF,  UKF,  and  UPF,  at  intervals  of  1  s,  10  s,  and  30  s  were  compared  with  the  control 
trajectory.  The  ranging  solution  updates  are  either  GPS  or  total  station  observations  and  all  three 
categories  of  inertial  sensors  were  tested.  Within  each  set  of  sensor/update  tests,  the  position 
errors  are  a  function  of  how  well  the  filter  estimates  the  systematic  errors  of  the  inertial  sensor  in 
the  particular  non-linear  environment. 

The  comparison  against  the  control  was  based  on  the  total-distance  errors,  computed  from  the 
three  coordinate  errors;  and,  their  standard  deviations  for  each  of  the  curved  and  straight  sections 
of  the  test  trajectory  (Figure  2)  were  averaged.  Figure  4  (for  1  s  updates  only)  shows  that,  as 
expected,  the  more  precise  updates  (total  station  vs.  GPS)  yield  more  accurate  free-inertial 
solutions  simply  because  the  integration  of  the  inertial  sensor  outputs  begins  with  a  smaller  error. 
Interestingly,  for  this  particular  simulation,  there  is  little  difference  between  the  high-end  and  the 
medium  quality  inertial  sensors  showing  that  either  one  will  offer  similar  interpolation  capability 
within  1  second.  Clearly,  the  EKF  does  not  perform  as  well  as  the  UT-based  filters,  particularly 
along  the  curved  segments  of  the  trajectory.  The  UPF  offers  only  slightly  better  performance 
than  the  UKF  in  these  tests  because  the  simulated  noise  processes  for  the  inertial  sensors  and  the 
observation  updates  are  Gaussian.  As  the  interval  before  the  next  update  increases,  the  inertial 
sensor  errors  accumulate  in  the  position  solution,  but  the  UKF  and  UPF  still  out-perform  the 
EKF,  as  shown  in  Figure  5  for  just  GPS  and  the  medium-grade  IMU.  We  note  that  the  position 
accuracy  of  the  filters  (no  additional  smoothing)  is  useful  only  for  screening  of  MEC  if  the  GPS 
outage  reaches  10  s. 
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Besides  outright  interruptions  in  the  ranging  solution,  it  may  also  be  degraded  during  several 
seconds  of  the  survey.  We  consider  only  differential  GPS  solutions  and  assume  that  such 
degradation  is  a  function  of  baseline  length  between  the  rover  and  the  fixed  base  station.  Using 
the  corresponding  increased  observation  noise,  as  shown  in  Table  2,  Figure  6  compares  the  EKF 
and  UKF  performances  with  respect  to  the  medium-grade  IMU,  free-inertial,  positioning 
accuracy  at  the  end  of  1  s  and  5  s  prior  to  degraded  GPS  updates.  Again,  we  find  improved 
results  for  the  UKF  over  the  EKF,  especially,  when  both  the  degradation  and  the  update  interval 
increase. 

The  superiority  of  the  UKF  may  be  realized  only  with  appropriate  tuning  of  the  scaling 
parameter,  a  (since  the  error  states  are  excited  by  Gaussian  noise,  the  optimal  value  J3  =  2  was 
used  for  these  tests).  The  results  of  Figure  4  and  5  were  obtained  with  a  =  0.15 ,  corresponding 
to  sigma  points  within  a  range  of  ±  lo  .  Expanding  this  range  did  not  improve  the  estimation, 
but  significantly  smaller  values  of  a  degrade  the  performance  of  the  UKF,  as  also  shown  in 
Figure  6,  which  includes  results  fora  =  0.0005 ,  a  =  0.001 ,  and  a  =  0.01 . 


Figure  4:  Standard  deviation  of  errors  of  different  IMU  grades  with  1  s  ranging  system  updates 
using  different  filters  (first  row:  curve  segments,  second  row:  straight  segments,  first  column: 
GPS,  second  column:  total  station)  (units:  cm). 
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Figure  5:  Standard  deviation  of  errors  for  different  filtering  methods  (first  row:  curve  segments, 
second  row:  straight  segments,  first  column:  10  s,  second  column:  30  s  GPS  updates  of  the 
medium-grade  IMU)  (units:  cm). 
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Figure  6:  Standard  deviation  of  errors  with  medium-grade  IMU  and  degraded  GPS  updates  and 
with  different  scales  in  the  unscented  transformation  (first  row:  curve  segments,  second  row: 
straight  segments,  first  column:  GPS  1  s  updates,  second  column:  GPS  5  s  updates). 


The  foregoing  tests  assumed  a  system  driven  by  Gaussian  noise.  We  may  expect  similar 
results  with  other  symmetrically  distributed  processes,  and  the  UKF  could  also  be  tuned  to  deal 
with  non-Gaussian,  symmetric  distributions  using  the  parameter,  |3  .  On  the  other  hand, 
dynamic  systems,  besides  being  generally  non-linear,  may  also  be  excited  by  non-symmetric 
processes  (Kushner,  1967).  Indeed,  Reddy  and  Herr  (2006)  investigated  the  skewness  of  IMU 
sensor  errors.  Julier  (1998)  and  Naveau  (2005)  also  proposed  modeling  the  process  noise  with 
asymmetric  probability  densities.  We  considered  both  symmetric  non-Gaussian  as  well  as 
asymmetric  distributions  to  determine  the  performance  of  the  particle  filter.  In  one  test,  sensor 
errors  were  generated  from  a  uniform  distribution  with  the  same  variance  as  in  the  Gaussian 
case.  As  seen  in  Table  3,  the  latitude  and  longitude  standard  deviations  are  slightly  lower  in 
value  with  the  particle  filter.  There  is  very  little  or  no  difference  with  respect  to  the  Gaussian 
case  (not  shown). 
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Table  3.  Statistics  of  errors  (all  segments)  using  nonlinear  filters  on  data  from  a  medium-grade 
IMU  with  uniformly  distributed  errors  and  with  1  Hz  GPS  updates  (units:  cm). 


coordinates 

UKF 

UPF 

avg. 

std. 

dev. 

avg. 

std. 

dev. 

0.1 

1.9 

0.1 

1.7 

A 

-0.02 

2.0 

-0.1 

1.9 

h 

0.6 

2.0 

0.5 

2.0 

To  test  the  efficacy  of  the  particle  filter  (PF),  we  considered  the  following  asymmetric 
probability  density  (APD)  (Komunjer,  2007). 
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where  0 <a  < 1 ,  A  >  0  and  oa  x  = - - - —  .  If  a  =  0.5 ,  then  the  density  is  symmetric.  If, 

ax+(l-a) 

in  addition,  A  =  2 ,  the  APD  is  the  Gaussian  density.  For  our  tests,  we  chose  A  =  2  and  two 
values  a  =0.25,0.75  (Figure  7),  one  for  the  gyro  noise,  the  other  for  the  accelerometer  noise. 
Figure  8  compares  the  statistics  of  the  position  error  for  the  update  rates  and  GPS  degradation 
used  in  Figure  6,  with  and  without  the  asymmetry  in  the  IMU  sensor  noises.  We  find  that  both 
the  UKF  and  EKF  performances  are  significantly  degraded  in  the  presence  of  IMU  noise 
asymmetry;  and,  moreover,  the  UKF  is  not  consistently  superior  to  the  EKF,  especially  in  the 
straight  sections  of  the  trajectory.  Even  in  the  curved  sections,  the  advantage  of  the  UKF  seen 
earlier  has  been  compromised  in  the  case  of  longer  ranging  gaps  and  a  degraded  GPS  solution. 
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pdf  for  gyro  noise  (a=0.25)  pdf  for  accelerometer  noise  (a=0.75) 


Figure  7:  Asymmetric  pdf’s  used  for  gyro  noise  (A,  =  2,  a  =0.25)  and  accelerometer  noise 
(A  =  2,  a  =0.75). 


On  the  other  hand,  the  unscented  particle  filter  (UPF)  applied  to  the  IMU  data  corrupted  by 
the  asymmetric  noise  exhibited  (Figure  9)  consistently  superior  performance,  both  for  the  curved 
and  the  straight  segments  of  the  trajectory.  This  is  an  interesting  result  because  the  UPF  still 
uses  the  unscented  Kalman  filter  before  the  re-sampling  of  the  particles;  yet,  the  results  are  better 
than  for  the  UKF,  alone. 

Figure  9  also  shows  that  increasing  the  number  of  particles  does  not  yield  significant 
improvements  in  the  UPF.  However,  by  incorporating  the  adaptive  algorithms  we  further 
enhance  the  filter,  at  least  for  the  longer  GPS  outage,  because  the  AUPF  compensates  for  the 
unfulfilled  assumptions  of  the  unscented  transformation,  specifically,  that  it  assumes  a  Gaussian 
(at  least  symmetric)  error  distribution.  The  overall  positioning  accuracy  deteriorates  as  the  GPS 
update  interval  increases  (e.g.,  5  s).  For  example,  the  standard  deviation  in  position  using  the 
UPF  with  asymmetric  error  distribution  is  worse  than  that  of  the  UKF  with  symmetrically 
distributed  errors.  However,  the  UPF  still  performs  better  than  the  other  filters  when  the  sensor 
noise  is  asymmetric. 
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Figure  8:  Standard  deviation  of  medium- grade  IMU  position  errors  with  asymmetric  sensor  error 
distributions,  for  baseline-degraded  GPS  updates  (first  row:  curve  segments,  second  row:  straight 
segments,  first  column:  1  s  updates,  second  column:  5  s  updates). 
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Figure  9:  Standard  deviations  of  medium-grade  IMU  position  errors  (first  row:  curved  segments, 
second  row:  straight  segments)  for  EKF,  UKF  and  UPF  with  1  s  GPS  updates  (first  column)  and 
for  EKF,  UKF,  UPF1  (200),  UPF2(400),  UPF3(600),  AUPF(200)  with  5  s  updates  (second 
column).  The  number  in  parenthesis  is  the  number  of  particles.  GPS  update  accuracy  is  for  the  1 
km  baseline  (Table  2)  and  the  IMU  sensor  noise  is  assumed  to  be  asymmetric  as  in  Figure  7 
(units:  cm). 
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4.  Performance  Using  IMU  Data 


As  part  of  this  project  we  purchased  two  medium  grade  IMUs,  the  HG1900  and  the  HG1700, 
from  Honeywell,  Inc.  The  HG1700  (Figure  10a)  is  a  light-weight  and  low-cost  ring-laser-based 
strapdown  IMU  designed  for  navigation  purposes  especially  for  missile  guidance  and  unmanned 
vehicles  (UAV).  The  HG1700  includes  three  miniature  GG1308  Ring  Laser  Gyroscopes  (RLGs) 
and  three  RBA-500  Resonating  Beam  Accelerometers.  The  Honeywell  HG1900  (Figure  10b)  is 
a  MEMS  (micro-electro-mechanical  system)  based  IMU  used  for  the  same  purposes  as  the 
HG1700,  but  slightly  less  accurate.  It  contains  MEMS  gyros  and  also  RBA500  accelerometers. 
A  comparison  of  the  technical  specifications  of  these  two  units  is  provided  by  the  manufacturer 
in  Figure  1 1 . 


Figure  10:  a)  HG1700  IMU,  b)  HG1900  IMU. 
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Figure  11:  Manufacturer  specifications  for  HG1700  and  HG1900. 


Prior  to  deploying  these  IMUs  in  a  geolocation  system  in  the  field,  we  tested  them  in  the 
laboratory  in  the  static  mode  for  noise  characteristics  and  on  a  cart  for  controlled  navigation. 
The  geolocation  filters  were  tested  both  in  the  laboratory  using  the  moving  cart  on  a  predefined 
trajectory  and  in  the  field  at  a  UXO  detection  test  site  in  collaboration  with  NRL  (Naval 
Research  Laboratory).  In  addition,  we  tested  a  new  self-calibration  procedure  for  the  IMUs 
locally  in  a  parking  lot  using  a  cart-based  system  with  GPS.  The  following  sections  provide  the 
details  of  these  tests. 
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4.1  Static  Laboratory  Tests 

4.1.1  Test  Description 

Data  from  the  HG1700  and  HG1900  units  were  collected  in  the  stationary  mode  (with 
presumably  constant  room  temperature)  for  the  maximum  time  allowed  by  the  data  acquisition 
software,  1000  seconds  at  100  Hz  (this  produces  100,000  data  records).  The  units  were  placed 
on  an  isolated  pillar  in  the  basement  of  Mendenhall  Laboratory  at  Ohio  State  University,  which 
is  meant  for  seismic  and  gravity  instrumentation.  Therefore,  the  units  are  unaffected  by  external 
vibrations  during  these  tests.  The  data  were  collected  using  a  run-box  and  software  purchased 
from  Honeywell. 

4.1.2  Results  and  Analysis 

The  raw  data  were  analyzed  by  subjecting  them  (unsmoothed)  to  a  Fourier  transform  to  compute 
their  power  spectral  density  (periodogram).  Figure  12  compares  the  psd’s  for  the  HG1700  and 
HG1900  accelerometers.  The  oscillations  evident  in  the  IMU  output  (not  shown)  are  clearly 
indicated  by  peaks  in  the  psd  at  0.02  Hz  and  possibly  0.04  Hz,  corresponding  to  50  s  and  25  s 
periods.  No  such  resonances  are  apparent  in  the  psd’s  of  the  HG1900  accelerometer  data. 
Moreover,  the  HG1900  psd’s  generally  have  much  lower  amplitude  at  the  lower  frequencies  (less 
than  1  Hz) 
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Figure  12:  PSD’s  of  the  HG1700  (top)  and  HG1900  (bottom)  accelerometers  in  the  static  mode. 

Similarly,  psd’s  were  generated  for  the  gyro  data  and  are  plotted  in  Figure  13.  The  psd’s  are 
similar,  although  those  for  the  HG1700  have  slightly  larger  amplitude  in  the  frequency  band  (0.1 
Hz  -  1  Hz)  and  significantly  greater  amplitude  for  higher  frequencies. 
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Figure  13:  PSD’s  of  the  HG1700  (top)  and  HG1900  (bottom)  gyros  in  the  static  mode. 


Table  4  shows  that  while  both  the  HG1700  and  HG1900  accelerometers  fall  within  specification 
(Figure  11),  the  HG1700  accelerometers  clearly  are  noisier.  Similarly,  for  the  gyros  we  find  that 
the  HG1700  gyros  are  slightly  noisier  (Table  5).  However,  the  gyro  biases  of  the  HG1700  (as 
determined  by  comparison  to  Earth’s  known  angular  rate)  are  somewhat  lower  (also  shown  in 
Table  5).  Both  the  biases  and  angular  random  walks  are  within  the  specifications  noted  in  Figure 
11. 
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Table  4:  Noise  levels  for  the  IMU  accelerometers  (static  mode). 


white  noise 

velocity  random  walk 

HG1700  accel 

4.8X10"4  m/sWHz 

2.9  xlO-2  (m/s)/Vhr 

HG1900  accel 

8.6xl0“5  m/sWHz 

5.2xl0“3  (m/s)/Vhr 

Table  5:  Noise  levels  for  the  IMU  gyros  (static  mode). 


bias 

angular  random  walk 

HG1700  gyro 

0.32  °/hr 

0.061  7x/hr 

HG1900  gyro 

1 

o 

C\ 

o 

i-i 

0.039  °/Vhr 

We  conclude  that,  overall,  the  HG1700  IMU  does  not  perform  as  well  in  the  static  mode  as  the 
HG1900.  Nevertheless,  when  using  these  units  to  navigate  a  cart  in  laboratory  tests  (and  in 
subsequent  field  tests),  we  consistently  obtained  slightly  better  positioning  results  with  the 
HG1700.  The  resonances  of  the  HG1700  accelerometers  at  0.02  Hz  and  0.04  Hz  have  an 
unknown  origin  (the  dither  of  the  gyros  required  to  prevent  the  lock-in  phenomenon  is  likely  not 
the  cause  since  the  dither  frequency  is  much  higher).  Working  briefly  with  the  technical  staff  at 
Honeywell,  we  were  not  able  to  ascertain  the  reason  for  the  degraded  noise  levels  of  the 
HG1700. 


4.2  Cart-Based  Laboratory  Tests 

We  tested  the  filter/smoothers  (AEKF  and  AUKF),  and  neural-network  aided  adaptive 
filter/smoothers  (NN-AEKS  and  NN-UKS)  in  the  laboratory  using  the  IMUs  mounted  on  a  cart. 
The  state  vector  included  three  position  errors,  three  velocity  errors  and  three  orientation  errors 
associated  with  the  navigation,  and  a  bias  and  a  scale  factor  error  for  each  of  the  three 
accelerometers  and  the  three  gyros.  For  the  UKF,  four  quaternions  were  employed  instead  of 
three  angles  for  the  orientation  error. 

In  addition  to  the  two  medium-grade  IMUs,  this  cart-based  geolocation  system  (CBGS)  also 
contained  the  H764G  navigation-grade  IMU  (similar  to  the  LN100  discussed  in  Section  3.1).  In 
addition,  the  cart  contains  the  IMU  data  collection  computer  hardware  and  a  physical  pointer 
used  to  identify  the  cart’s  passing  or  occupation  of  a  control  point  (Figure  14).  The  pointer 
served  a  function  similar  to  an  external  position  observation  (such  as  from  GPS).  The  error 
specifications  of  the  three  IMUs  as  provided  by  the  manufacturer  are  described  in  Table  6. 
These  specified  error  parameter  values  are  used  to  generate  the  initial  covariance  matrix  of  the 
system  noise  ( Q ,  see  equation  (4)).  The  accelerometer  bias  of  the  H764G  is  about  50  times 
smaller  than  for  the  medium-grade  units  and  its  gyro  biases  are  similarly  much  smaller. 
Although  the  static  tests  indicated  somewhat  poorer  noise  characteristics  for  the  HG1700,  we 
assigned  the  manufacturer-specified  variances  to  the  Q-matrix  and  initial  P- matrix  (for  the 
biases)  (see  equation  (4)). 
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Table  6:  The  error  specification  of  three  IMUs 


H764G 

HG1700 

HG1900 

Accel 

Error 

Bias 

20  Mg 

1  mg 

1  mg 

Scale  Factor 

40  ppm 

300  ppm 

300  ppm 

Random  Walk 

0.003  (m/s)/Vhr 

0.09  (m/s)/Vhr 

0.09  (m/s)/Vhr 

Gyro 

Error 

Bias 

0.01  deg/hr 

1  deg/hr 

<7  deg/hr 

Scale  Factor 

1  ppm 

150  ppm 

150  ppm 

Random  Walk 

0.001  deg/Vhr 

0.125  deg/Vhr 

0.09  deg/Vhr 

Figure  14:  Cart  Based  Geolocation  System,  Front  (A:  HG1700  and  HG1900  with  run-box,  B 
H764G,  C:  Pin  Point  Indicator  with  Mark) 


4.2.1  Test  Description 

In  our  laboratory  tests,  we  used  simulated  position  updates  (a  manual  coordinate  registration 
system)  instead  of  integrating  the  IMU  with  a  ranging  system  such  as  GPS.  The  cart  was  pushed 
along  a  trajectory  with  pre-defined  waypoints  that  have  known  coordinates.  Whenever  the  cart 
passed  a  particular  waypoint,  the  event  was  recorded  with  a  timing  mark  from  the  computer 
clock  that  also  time-tagged  the  collection  of  data  from  the  IMU.  The  imperfect  recording  of  the 
passage  of  the  cart  pointer  (Figure  14C)  over  the  ground  marker  could  be  considered  an  error  in 
the  control  point  coordinates,  although  it  is  a  personal  error  and  not  as  large  in  magnitude  as  a 
kinematic  GPS  position  error.  The  magnitude  of  this  error  is  estimated  to  be  less  than  1  cm  per 
control  point;  however,  the  variance  of  the  measurement  error  (diagonal  matrix,  R  ,  see  equation 

(6))  was  conservatively  set  to  (lcm)2 .  The  test  trajectory  had  four  straight  lines  and  four  curved 
sections  (Figure  15).  Twenty-eight  (28)  ground  marks  were  used  as  “measured”  control  points 
for  the  filtering/smoothing.  The  distance  between  points  is  12  inches  and  the  CBGS  needed 
about  4  seconds  to  move  from  one  control  point  to  another.  Therefore,  the  speed  of  the  cart  is 
about  0.076  m/s. 
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Figure  15:  Laboratory  test  trajectory. 


We  tested  the  integrated  system  by  comparing  positioning  solutions  at  control  points  not  used  in 
the  integration.  In  the  first  case,  every  other  point  was  used  to  update  the  system  and  the 
filtered/smoothed  position  estimate  was  compared  to  the  skipped  control.  The  total  duration  of 
inertial  positioning  between  control  updates  is  about  8  seconds  in  this  situation.  In  the  second 
case,  every  third  point  served  as  update  and  the  estimated  positions  were  compared  to  the 
coordinates  of  the  two  skipped  points.  The  total  simulated  GPS  outage  in  this  case  is  about  12 
seconds.  The  (neural-network  aided)  filtered/smoothed  position  errors  were  analyzed  along  the 
straight  and  curved  sections,  respectively,  in  terms  of  the  standard  deviation  of  the  total  3-D 
position  error  computed  from  5  separate  tests  along  the  same  trajectory.  In  the  first  case,  there 
are  8  (6)  comparison  points  in  the  straight  (curved)  segments,  hence  40  (30)  comparisons 
contribute  to  the  standard  deviation  of  the  error.  In  the  second  case,  the  corresponding  number 
of  comparisons  is  45  for  both  straight  and  curved  segments. 


4.2.2  Results  and  Analysis 

Figures  16  and  17  show  the  standard  deviations  of  position  errors  of  the  CBGS  according  to  the 
different  IMUs,  the  various  filtering/smoothing  methods,  and  the  interval  between  control 
updates.  Figures  18  and  19  show  the  corresponding  standard  deviations  when  the  neural- 
network  (NN)  aided  adaptive  filter/smoothers  were  applied.  As  anticipated,  the  navigation-grade 
IMU  (HG764G)  performs  best.  Also,  the  non-linear  based  smoothing  method  (AUKS)  yields 
better  results  than  the  EKF-based  smoothing  method,  especially  in  the  turning  segments.  The 
HG1900  performs  only  slightly  worse  than  the  HG1700.  The  relative  difference  in  performance 
between  the  tactical  or  medium  grade  IMUs  (HG1700  and  HG1900)  and  the  navigation-grade 
IMU  in  these  tests  decreased  significantly  using  the  non-linear  filtering  techniques. 
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Figure  16:  Position  Errors  of  CBGS  according  to  different  IMUs  and  filtering/smoothing 
methods  (control  updates  every  2  points). 
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Figure  17:  Position  errors  of  CBGS  according  to  different  IMUs  and  filtering/smoothing  methods 
(control  updates  every  third  point). 


Similar  to  the  two-point  update  case,  in  every  three-point  update  case  the  H764G  achieved 
the  best  position  accuracy;  and,  the  nonlinear,  adaptive  filter/smoothing  techniques  demonstrated 
better  performance  than  the  AEKS.  However,  compared  to  the  first  case,  there  was  less 
difference  with  both  method  (AEKS  and  AUKS)  between  the  straight  and  turning  segments 
because  the  accumulation  of  IMU  errors  overwhelmed  the  superiority  of  the  nonlinear  filter  in  a 
dynamic  environment.  On  the  other  hand,  the  non-linear  filter  still  achieves  better  estimates  of 
states  in  an  essentially  non-linear  system. 

The  overall  position  accuracy  was  improved  by  including  the  neural  network  aiding,  and  the 
disparity  in  performance  between  straight  and  turning  sections  was  reduced.  It  is  noted  that  the 
position  error  of  the  NN-AEKS  decreased  dramatically  compared  to  the  AEKS  in  both  straight 
and  curved  section.  The  H764G  with  NN-aided,  adaptive  filtering/smoothing  can  achieve  the 
area-mapping  position  requirement  (5  cm)  along  straight  sections  with  8  seconds  between 
updates.  Similar  results  were  obtained  with  the  longer  interval  between  updates  (every  third 
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control  point).  The  position  error  (st.  dev.)  of  the  H764G  is  less  than  10  cm  in  this  case,  but  in 
particular,  the  HG1700  and  HG1900  yielded  only  slightly  worse  results  (l~3cm  more  in  standard 
deviation).  With  the  increased  interval  between  control  points  the  NN-aided,  adaptive 
filter/smoother  was  able  to  maintain  the  position  accuracy  better  than  the  unaided  adaptive 
filter/smoother. 
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Figure  18:  Position  Error  of  CBGS  according  to  different  IMUs  and  NN  aided 
filtering/smoothing  method  (control  updates  every  2  points). 


Figure  19:  Position  Error  of  CBGS  according  to  different  IMUs  and  NN  aided 
filtering/smoothing  method  (control  updates  every  3  points). 


4.3  Hand-Held  System  Tests 

The  Handheld  Geolocation  System  (HGS)  developed  for  our  tests  had  a  tactical-grade  IMU 
(HG1900)  and  was  integrated  with  an  automatic  target  tracking  system  that  updates  its  position 
(Figure  20)  at  a  specified  rate.  The  manufacturer  error  specifications  of  the  HG1900  are 
described  in  Figure  1 1  and  Table  6.  A  red  circle  marked  on  the  front  head  of  the  HGS  was  used 
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as  a  target  that  could  be  tracked  for  accurate  positioning  by  using  real-time  color  imaging.  This 
external  tracked  position  served  the  same  purpose  as  would  GPS  in  the  field. 

In  addition,  we  developed  improved  software  to  collect  the  data  from  the  IMU,  thus  replacing 
the  Honeywell  software  (that  was  an  MS-DOS  based  data  acquisition  tool)  that  limited  the 
number  of  collected  data  to  100000  records.  This  also  enabled  us  to  collect  data  using  a  laptop 
(replacing  the  desktop  computer  used  with  the  cart-based  system).  The  run  box  obtained  from 
Honeywell  still  serves  as  a  power  and  data  cable  terminal,  but  otherwise  has  no  function  in  the 
data  collection. 


Figure  20.  The  hardware  system  for  locating  geolocation  tests  (A:  Handheld  Geolocation 
System,  B:  Laptop,  C:  IMU  Run-box,  D:  Web  Cam,  E:  Target  (red  dot)  with  black  box,  F: 
HG1900,  G:  Target  tracking  software,  and  H:  PCMCIA  converter). 


In  order  to  obtain  accurate  positions  of  the  HGS,  an  automatic  target  position  extraction 
system  was  implemented  using  real-time  color  imaging  with  a  webcam  (at  5Hz  update  rate). 
First,  the  color  image  is  converted  to  a  binary  image  which  has  only  two  values,  0  or  255,  for 
each  R,  G,  and  B  band  (only  the  R  band  was  employed  in  this  test).  Second,  the  target  object  in 
the  binary  image  is  extracted  from  the  boundary  noise  using  geometric  circularity  and  knowledge 
of  its  size.  Finally,  the  screen  coordinates  of  (the  center  of  the  circle)  are  converted  to  real  world 
coordinates. 

Using  our  developed  target  position  tracking  software  the  transformed  binary  image  is 
captured  and  displayed  in  the  right-top  section  of  the  laptop  display  (Figure  21)  with  a  pre¬ 
defined  but  adjustable  threshold  value  determined  with  a  slide  bar.  The  minimum  size  of  the 
target  and  the  circularity  to  differentiate  the  target  from  other  objects  are  defined  as  100  and  5, 
respectively.  Since  the  actual  map  size  (128cm*960cm)  will  be  used  to  convert  the  target 
position  in  the  window  coordinate  system  (640x480  [pixel])  to  the  position  in  the  real-world 
system,  the  actual  dimension  of  one  pixel  is  0.2cm  by  0.2cm.  The  position  of  image  plane  of  the 
camera  (leveled  using  bubble  levels)  and  the  test  area  (laboratory  floor)  are  assumed  parallel, 
with  the  camera  located  1  m  above  the  center  of  the  test  area.  An  error  of  less  than  0.5  degree 
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due  to  the  misalignment  between  the  camera  and  the  test  area  implies  a  final  position  error  that  is 
less  than  0.87cm  (  =  0.5°  xlm ,  or  bout  4  pixels). 


4.3.1  Test  Description 

The  HGS  was  tested  over  a  small  area  with  its  position  tracked  by  the  web  cam  and  tracking 
software,  as  described  in  the  previous  section.  The  RBPF  and  the  previously  developed 
filters/smoothers  (KF/KS  and  UPF/UPS)  were  evaluated  for  typical  local  trajectories  of  a  hand¬ 
held  UXO  detection  system.  Without  loss  in  generality,  the  loosely  coupled  INS/GPS 
integration  scheme  based  on  the  decentralized  filter  architecture  was  employed.  The  state  vector 
for  the  IMU/imaging  system  comprised  19  states:  two  horizontal  position  errors  (dxN,dxE)  in  the 
navigation  frame;  two  velocity  errors  (8vw,8v£),  three  orientation  errors  in  a  local  north-east- 
down  frame  (\|/W,\)/£,\|/D);  and,  biases  ( bA  and  bG )  and  scale  factor  errors  (k4  and  kg)  for  the 
accelerometers  and  gyros. 

The  state  vector  can  be  divided  using  the  Rao-Blackwellization  described  in  Section  2.2.3 
into  two  parts  according  to  Nordlund  (2002)  and  Hektor  (2007).  They  argued,  based  on 
simulations  and  airborne  test  data,  that  only  the  position  states  are  highly-nonlinear  and  all  other 
states  can  be  assumed  as  linear  without  significant  loss  in  position  accuracy.  However,  this  state 
division  is  incorrect  in  our  case  because  not  only  the  position  error  states  but  also  the  orientation 
and  velocity  error  states  experience  high  nonlinear  dynamics  in  a  ground-based  UXO  detection 
system  such  as  the  Hand-held  Geolocation  system.  Therefore,  the  division  of  states  should  be 
separated  into  the  navigation  related  states  (position  errors,  velocity  errors  and  orientation  errors) 
and  all  other  states  (bias  and  scale  factor  error  of  the  gyros  and  accelerometers): 

(xj.)  —  [8xw  8x£  8 vN  8v£  \\t n  V|/£  t|/D]  ,  (37) 

(4)T  =[K  K  K  b x  k<4  kCi  K<i  k4  k4  k.J1',  os) 

The  typical  dynamics  of  a  hand-held  UXO  detection  platform  can  be  classified  into  four 
categories;  linear,  curved,  sweep,  and  swing  (Bell  and  Collins,  2007).  However,  since  the  linear 
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and  curved  sections  are  included  in  the  sweep  and  swing  dynamics,  only  two  test  trajectories 
(sweep  and  swing)  were  considered.  The  position  accuracy  of  the  HGS  system  was  tested  along 
five  sweep  and  five  swing  trajectories.  The  sweep  trajectory  had  six  straight  lines  and  five 
curved  sections.  The  swing  trajectory  had  five  straight  lines  and  four  curved  sections.  The  total 
distance  of  sweep  (swing)  of  the  trajectory  is  about  7.2  (5.6)  m  and  the  HGS  took  about  22  (14) 
seconds  to  complete  total  trajectory.  Therefore,  the  speed  of  the  HGS  was  about  0.33  (0.4)  m/s, 
respectively.  Figure  22  shows  the  trajectories  of  the  first  sweep  and  swing  of  five  separate  tests. 


X  (cm)  X  (cm) 

Figure  22.  The  first  sweep  and  swing  trajectory  from  five  tests  of  the  designed  handheld  UXO 
Geolocation  platform. 


4.3.2  Test  Results  and  Analysis 

Various  update  intervals  of  the  imaged  target  positions  were  implemented  to  analyze  the 
filtered/smoothed  position  estimates  with  respect  to  the  requirements  of  position  accuracy  for 
MEC  geolocation  and  characterization.  In  each  test  case,  5Hz,  2.5Hz,  1Hz,  0.5Hz,  and  0.25Hz 
update  rates  were  employed  in  the  integration.  The  accuracy  of  the  HGS  was  tested  by 
comparing  the  estimated  position  to  the  target  positions  (available  at  5Hz)  which  were  not  used 
as  updates  in  the  filtering  process.  For  example,  if  the  update  rate  of  the  filter  is  1Hz,  every  fifth 
point  of  the  automatically  tracked  target  position  is  used  as  external  control  (the  accuracy  of 
measurement  is  assumed  as  2~4  pixel  (0.4cm)  ~(0.8cm)  in  north  or  east  direction)  in  the  filter 
and  the  other  four  points  are  compared  to  the  estimated  points  of  the  filter/smoother.  The 
standard  deviations  of  the  total  2-D  errors  (from  the  EKS,  UPS,  and  RBPS)  were  computed  for 
the  straight  and  curved  sections.  The  curved  section  is  defined  by  a  pre-defined  threshold  of 
absolute  value  of  heading  change  angle  ( >  20°). 

Sweep  Test:  Figure  23  shows  the  average  standard  deviations  of  the  error  for  each  of  the 
filters  implemented  and  for  each  update  rate  of  the  control  points.  The  averages  represent  the 
results  of  the  five  separate  sweep  tests.  In  the  straight  section,  the  UPS  and  RBPS  can  achieve 
the  discrimination  level  of  position  accuracy  (better  than  2cm  standard  deviation)  for  update 
rates  ranging  from  5Hz  to  a  0.5Hz,  and  they  can  achieve  the  area  mapping  level  of  accuracy 
(better  than  5cm  for  standard  deviations)  at  the  0.25Hz  update  rate.  The  EKS  yields  comparable 
(slightly  worse)  results  at  0.5Hz  update  rate  and  significantly  worse  results  when  the  update  rate 
is  0.25Hz.  In  the  curved  section,  every  filter  approached  the  discrimination  level  of  accuracy 
(less  than  2cm  std.  dev.)  up  to  1Hz  update  rate  and  achieved  5cm  (std.  dev.)  accuracy  at  0.5Hz. 
The  position  errors  of  the  RBPF  are  comparable  with  the  UPS  up  to  0.5Hz  but  slightly  lower 
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than  all  other  smoothers  at  0.25Hz.  Therefore,  the  RBPS  is  expected  to  provide  superior 
performance  when  the  external  imaging  solution  is  blocked  or  degraded  for  longer  than  4 
seconds. 
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Figure  23.  The  average  standard  deviations  of  position  errors  according  to  various  smoothing 
methods  (unit:  cm)  for  the  sweep  tests. 


Swing  Test :  Overall,  the  results  of  the  swing  test  were  worse  than  those  from  the  sweep  test, 
especially  when  the  update  rate  is  less  than  1Hz.  Similar  to  the  sweep  test,  the  non-linear 
filter/smoothing  techniques  demonstrated  better  performance  than  the  EKS.  However,  compared 
to  the  sweep  test,  there  was  less  of  a  difference  between  straight  and  curved  sections  in  all 
estimation  methods.  Figure  24  shows  the  average  standard  deviations  of  position  errors 
according  to  the  various  smoothing  methods  and  the  update  interval  between  control  points.  The 
position  accuracy  of  the  swing  test  was  degraded  compared  to  the  sweep  test,  especially  when 
the  update  rate  is  less  than  1Hz  in  both  straight  and  curved  sections.  However,  it  is  noted  that  the 
position  error  of  the  UPS  and  RBPS  is  significantly  smaller  at  0.25Hz  when  compared  to  the 
EKS  in  the  straight  and  curved  sections. 

In  the  straight  section,  the  UPS  and  RBPS  attained  the  discrimination  level  of  position 
accuracy  up  to  the  1Hz  update  rate  and  achieved  the  area  mapping  level  of  accuracy  at  0.5Hz 
update  rate  (compared  to  the  sweep  test  where  the  update  rates  for  the  discrimination  and  area 
mapping  level  of  accuracy  were  0.5Hz  and  0.25Hz,  respectively).  In  the  curved  section,  similar 
to  the  straight  section,  the  UPS  and  RBPS  achieved  the  discrimination  level  up  to  1Hz  and  the 
area  mapping  level  up  to  0.5Hz. 

In  this  swing  test,  the  RBPS  performs  best  at  all  update  rates.  Therefore,  although  the  swing 
operation  can  obtain  position  data  for  the  UXO  detection  in  shorter  time  duration  and  therefore 
with  fewer  control  points,  the  sweep  operation  is  preferred  because  the  position  errors  of  the 
sweep  test  were  smaller  than  that  of  the  swing  test  due  to  the  linear  dynamics  of  the  sweep 
motion. 
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Figure  24.  The  average  standard  deviations  of  position  errors  according  to  various  smoothing 
methods  (unit:  cm)  for  the  swing  tests. 


The  UPS  needs  at  least  200  particles  to  yield  optimal  position  accuracy  in  this  test.  However, 
the  RBPS  utilized  only  20  samples  (particles)  for  the  nonlinear  (particle)  filter  part.  Although 
there  is  still  room  for  some  position  accuracy  improvement  by  increasing  the  number  of 
particles,  it  will  not  yield  significant  improvements  (Lee  and  Jekeli,  2009).  Therefore,  the  RBPS 
can  produce  (slightly)  better  or  comparable  results  compared  to  the  UPS  with  only  10%  of  the 
number  of  particles  used  by  the  UPS. 


4.4  Field  Tests  at  NRL  Test  Facility 

The  filtering  methods  described  in  Section  2  were  evaluated  using  IMU/GPS  data  from  field  tests 
conducted  at  NRL’s  UXO  test  site  located  at  the  Army  Research  Laboratory  Blossom  Point 
Facility  in  Maryland  with  the  Multi-sensor  Towed  Array  Detection  System  (MTADS,  Figure 
25a)  and  a  cart-based  system  (Figure  26).  At  the  center  of  MTADS  we  mounted  OSU’s  dual- 
IMU/GPS  system  (Figure  25c):  two  tactical-grade  IMUs  (HG1700  and  HG1900)  in  the  IMU  box 
(Figure  25b),  one  Topcon  geodetic  GPS  receiver  (GB-1000),  and  a  laptop  computer  for  data 
collection  and  communication  with  the  sensors.  The  cart-based  system  has  only  OSU’s 
IMU/GPS  system  and  one  large  12V  battery  for  power  (Figure  26). 
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Figure  25.  (a)  The  NRL  Multi-sensor  Towed  Array  Detection  System  (MTADS),  (b)  IMU  box, 
(c)  OSU’s  dual-IMU/GPS  system. 


Figure  26.  The  NRL  cart-based  system  with  OSU’s  dual-IMU/GPS  geolocation  system. 


4.4.1  Test  Description 

Eight  field  tests  were  performed  using  the  two  platforms  and  various  GPS  satellite 
configurations.  Among  these  tests,  GPS  malfunctioned  for  one  complete  test,  and  three  tests 
were  designed  for  deliberate  GPS  degradation  (near  a  wooded  area).  Thus,  we  considered  only 
the  four  Tests  2,  5,  6,  and  7  in  this  analysis.  They  were  divided  into  two  test  scenarios  (Scenario 
1\  Tests  2  and  5,  Scenario  2:  Tests  6  and  7)  according  to  the  number  of  visible  GPS  satellites. 
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The  GPS  PDOP  (Position  Dilution  of  Precision)  for  Scenario  1  was  lower  than  that  for  Scenario 

2. 

The  test  range  is  roughly  84m  long  by  24m  wide.  The  vehicle-towed  trailer  and  cart-based 
system  were  moved  along  trajectories  that  are  typical  of  a  UXO  survey.  At  the  end  of  each 
trajectory  segment  there  are  large  turns  to  align  the  platform  to  the  next  straight  segment  of  the 
trajectory. 

Figure  27  shows  the  GPS  trajectories  of  the  vehicle-towed  trailer  and  the  cart-based  system 
for  Scenario  1  (Test  2,  MTADS;  and  Test  5,  Cart)  and  Scenario  2  (Tests  6,  Cart,  and  7, 
MTADS).  The  total  distance  of  Test  2  (Test  5)  is  about  2,016m  (398m)  and  the  speed  of  the 
MTADS  (Cart)  is  about  1.85  m/s  (0.95m/s).  The  total  distance  of  Test  6  (Test  7)  is  about  352m 
(514m)  and  the  speed  of  the  Cart  (MTDAS)  is  about  1.0  m/s  (1.9  m/s). 
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Figure  27.  The  GPS  trajectory  of  the  vehicle-towed  trailer  and  cart-based  System  (Tests  2  and  5 
are  Scenario  1  and  Tests  6  and  7  are  Scenario  2). 


The  kinematic  GPS  analysis  package  "KinTOOLS"  was  used  to  process  1Hz  GPS  data 
collected  by  the  Topcon  GB-1000  receiver.  Each  trajectory  using  the  individual  inertial  sensors 
was  estimated  by  the  EKF  or  the  UKF  and  then  the  wave  correlation  filter  was  applied  to  the 
common  solutions  with  the  threshold  value  of  0.5  (equation  (30)).  The  remaining  bias  and  trend 
errors  are  further  removed  by  the  end-matching  method.  Figure  28  illustrates  the  general 
IMU/GPS  data  processing  and  analysis  procedure. 
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4.4.2  Results  and  Analysis 

Scenario  1 :  The  filtered  free- inertial  (IMU-only)  positions  of  the  vehicle-towed  trailer  and  cart- 
based  system  were  computed  along  the  straight  sections  within  2,  4,  and  6  second  intervals 
representing  artificial  periods  of  GPS  unavailability.  Figure  29  illustrates  parts  of  the  estimated 
positions  from  the  trajectory  of  Test  2. 

As  the  interval  between  GPS  updates  increases  from  2  to  6  seconds,  the  IMU-estimated 
position  (HG1700)  using  the  EKF  with  wave  correlation  filter  and  end- matching  is  demonstrably 
better  than  the  position  using  just  the  EKF  and  the  EKF  with  WCF.  Similar  improvement  was 
obtained  with  the  WCF  applied  to  the  UKF  solutions. 
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2  sec  Update 


Figure  29.  The  estimated  vehicle  trajectory  using  GPS,  the  HG1700  w/  UKF  the  WCF,  and  end¬ 
matching,  for  different  update  rates. 


Figure  30  shows  the  standard  deviations  of  position  errors  of  the  three  coordinates  based  on 
all  straight  sections  of  the  trajectory)  of  Tests  2  and  5  according  to  the  UKF,  the  WCF  estimates 
based  on  the  UKF  solutions,  and  the  estimates  obtained  by  applying  the  UKF,  the  WCF,  and 
end-matching.  We  show  only  the  UKF-only  estimates  for  the  HG1700  because  it  yields  better 
position  accuracy  than  the  HG1900. 

Test  5  (cart-based  system)  yielded  overall  slightly  better  performances  than  Test  2  (vehicle- 
towed  trailer)  because  the  cart-based  system  experienced  lower  dynamics  and  slower  speed  than 
the  vehicle-towed  trailer.  The  WCF  with  end-matching  improved  the  UKF  solutions  (control 
updates  every  2  points)  of  Test  2  by  about  46%,  compared  to  an  11%  decrease  in  the  standard 
deviation  of  errors  without  the  end-matching. 

In  the  4  and  6  second  GPS  outage  cases,  the  position  errors  (standard  deviations)  decreased 
64%  and  76%  with  respect  to  the  UKF-only  errors,  and  55%  and  70%  relative  to  the  UKF  with 
WCF.  In  Test  5,  the  position  errors  with  the  WCF  and  end-matching  decreased  about  52  %  ~  69 
%  with  respect  to  the  UKF-only  solutions  and  46  %  ~  64  %  with  respect  to  the  UKF  plus  WCF 
solutions,  considering  all  2,  4,  and  6  second  GPS  outages. 
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Figure  30.  The  standard  deviation  of  position  error  in  the  straight  section  according  to  the  UKF, 
the  UKF  with  WCF,  and  the  UKF  with  WCF  and  end-matching  based  on  the  UKF  solutions 
(left:  Test  2,  right:  Test  5). 


Scenario  2:  In  Scenario  2  (Tests  6  and  7),  the  number  of  GPS  satellite  increased  from  4~5  to 
7-8,  and  the  PDOP  correspondingly  decreased  from  4.3  (Test  6)  to  2.6  (Test  7).  Similar  to 
Scenario  1,  the  WCF  and  end- matching  improved  the  performances  of  the  UKF  in  all  straight 
sections  (Figure  31).  In  Test  6  (Test  7),  the  position  errors  with  the  WCF  decreased  about  52  % 
-  69  %  (50  %  -  73  %)  with  respect  to  the  UKF-only  solutions,  and  44  %  -  63  %  (43  %  -  67  %) 
with  respect  to  the  UKF  with  WCF  solutions,  again  considering  all  GPS  outages. 


Figure  31.  The  standard  deviation  of  position  error  in  the  straight  section  according  to  the  UKF, 
the  UKF  with  the  wave  correlation  filter,  and  the  UKF  with  the  wave  correlation  filter  and  end¬ 
matching  based  on  the  UKF  solutions  (left:  Test  6,  right:  Test  7). 


4.5  Self-Calibration  Tests 

IMU  Calibration  is  the  process  of  determining  the  errors  in  the  outputs  of  the  gyros  and 
accelerometers  and  is  essential  in  order  to  increase  the  accuracy  of  the  navigation  solution.  In  the 
calibration  process,  the  deterministic  errors  (biases,  scale  factor  errors  and  sensor  non-orthogonality 
errors)  are  estimated  by  comparing  the  instrument  outputs  with  known  data,  and  the  random  errors 
(noise)  are  minimized  by  filtering  methods.  Various  kinds  of  calibration  methods  for  determining  the 
IMU  errors  have  been  developed.  The  most  common  methods  use  precise  laboratory  instruments, 
but  method  based  on  field  data  can  also  lead  to  reasonably  accurate  calibration,  particularly  methods 
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that  compare  the  outputs  in  the  static  mode  to  natural  signals  generated  by  the  Earth:  its  gravitational 
attraction  and  rotation  rate. 

We  tested  the  multi-position  calibration  method  developed  by  Shin  (2001)  on  the  two 
medium  level  IMUs:  HG1700  and  HG1900.  Because  the  Earth’s  rotation  rate  is  relatively  weak, 
some  of  the  systematic  errors  in  the  gyros  could  not  be  estimated;  nevertheless,  we  assumed  a 
full  model  of  bias,  scale  factor  error,  and  non-orthogonality  for  all  six  sensors  of  the  IMU.  The 
equations  for  the  accelerometer  and  gyroscope  measurements,  ya  and  y0) ,  respectively,  are 
therefore  given  by 


v  =a+b+Sa  +lf/ aa  +  wa , 
v  =?  +  h  -f\  ?  4-lJ/  ?  +w 

■T«>  •  ^  Jfi>'  T  "co  ’ 


(39) 

(40) 


where  a  and  ?  are  the  true  acceleration  and  angular  rate  vectors;  ba ,  bm ,  are  the  respective  bias 
vectors;  Sa,  Sm,  are  the  scale  factor  error  matrices;  ( ,  *T  0)  are  the  non-orthogonality  error 
matrices;  and  wa ,  wm  are  the  noise  vectors. 

For  field  operations,  where  the  IMU  is  assumed  capable  of  being  oriented  in  arbitrary 
directions,  we  tested  the  multi-position  calibration  method  developed  by  Shin  and  El-Sheimy 
(2002);  see  also  Syed  et  al.  (2007).  By  recording  the  IMU  data  with  different  attitudes  in  the 
static  mode,  the  deterministic  sensor  errors  are  estimated  to  the  extent  possible  by  the  natural 
signals.  For  a  given  attitude,  the  calibration  model  in  the  sensor  frame  for  the  accelerometers  is 
given  by 


fa  =ax+a+a  -  g  =0, 


(41) 


where  g  is  the  gravity;  and,  the  true  accelerations  are,  from  equation  (39), 


y  —  b 
a  =3^ - 

x  h  ,  ax 

l  +  Sa, 

...  ,  y ay -Ky 

ay=aWa,yz+— - ~~W- 


l  +  s„ 


ay 


(42) 


y  —b 

a  =  —a  \l/  +  a  \l/  H — — - —  —  w 

z  XT  a,zy  jt  a,zx  X ^  ° 


where  the  notation  is  self-evident.  Similarly,  for  the  gyroscopes  the  model  in  the  sensor  frame  is 

(43) 


/co  =°>*2  +  cov2  +  coz2-(D£  2  =0, 


where  (tiF  is  Earth’s  rotation  rate.  From  equation  (40),  we  find  analogously 
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These  models  were  applied  to  the  outputs  for  different  orientations  of  the  sensor  and  the  error 
parameters  (biases,  scale  factor  errors,  non-orthogonality  errors)  were  solved  by  non-linear  least- 
squares  adjustment.  For  details,  see  Hayal  (2010). 


4.5.1  Simulation  Tests 

The  model  described  in  the  previous  section  was  tested  using  simulated  measurements  of  an 
IMU  such  as  the  HG1700  with  noise  simulated  according  to  parameters  given  in  Table  6.  The 
data  were  generated  by  mathematically  orienting  the  sensors  relative  to  the  Earth  signals  with  26 
different  attitudes,  where,  imaging  the  IMU  as  a  cube,  each  of  the  6  faces,  8  corners,  and  12 
edges  in  turn  was  pointed  down.  Referring  to  Hayal  (2010)  for  the  detailed  analyses,  these 
simulations  showed  that  with  these  measurements,  all  9  accelerometer  errors  (biases,  scale  factor 
errors,  and  non-orthogonality  errors)  can  be  estimated  accurately  (limited  only  by  noise).  Also, 
the  gyro  biases  could  be  estimated,  but  it  is  possible  to  estimate  accurate  gyroscope  scale  factor 
errors  and  non-orthogonality  errors  only  if  we  increase  the  reference  signal  strength  (e.g.,  several 
orders  of  magnitude  greater  than  Earth  rotation)  or  correspondingly  decrease  the  noise  in  the 
data.  In  this  case  the  bias-only  model  may  be  more  appropriate.  However,  our  field  tests  (see 
the  following  section)  indicated  that  this  may  not  necessarily  lead  to  better  results. 


4.5.2  Field  Tests 

To  analyze  the  multi-position  calibration  method  in  the  field  and  how  it  affects  the  positioning 
performance  of  an  integrated  IMU/GPS  geolocation  system,  we  instrumented  a  cart  with  the 
Honeywell  HG1700  and  HG1900  IMUs,  as  well  as  a  GPS  receiver  (Figure  32),  and  conducted 
the  appropriate  tests  in  a  campus  parking  lot  of  the  Ohio  State  University.  Initially  20  minutes  of 
warming  up  time  preceded  the  measurements.  As  done  in  the  simulations,  26  different  5-minute 
long  static  IMU  attitude  measurements  were  collected  for  estimation  of  the  IMU  bias,  scale 
factor  error  and  non-orthogonality  errors.  Then,  the  IMU  box  was  installed  on  the  cart  together 
with  the  GPS  receiver  and  the  remaining  instruments.  The  GPS  measurement  and  positioning 
intervals  were  set  equal  to  0.1  second.  Moreover,  another  GPS  receiver,  the  base  station,  was 
installed  at  a  local  control  point  and  its  measurement  interval  was  set  to  0. 1  second.  After  that, 
the  INS  cart  was  pulled  by  hand  about  25  minutes  with  a  walking  speed,  and  a  sweep  like 
trajectory  (see  Figure  33)  was  followed  among  the  parking  space  lines.  The  trajectory  included 
25  straight  sections  and  25  very  short  and  sharp  curved  sections. 
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Figure  32:  Cart-Based  Geolocation  System  with  1)  IMU  box  containing  the  HG1900  and 
HG1700;  2)  120  VAC  inverter;  3)  12  VDC  battery;  4)  hand-pulled  cart;  5)  Trimble  NERTS  GPS 
receiver;  6)  Trimble  Zephyr  Geodetic  II  GPS  antenna;  7)  laptop  computer,  including  PCMCI 
card  connecting  the  IMUs  and  the  decoding  software. 
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Figure  33:  GPS  trajectory  of  field  IMU  calibration  tests. 


Various  combinations  of  attitudes  (of  the  26  available)  and  parameters  (of  the  9  for  each 
accelerometer  or  gyro  triad)  were  tested  to  determine  which  scenario  would  yield  the  best 
estimates  in  the  least-squares  adjustment.  Using  just  the  first  10  attitudes  gave  the  smallest 
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internal  standard  deviation  in  the  adjustment  for  the  gyro  biases  and  accelerometer  biases,  scale 
factor  errors  and  non-orthogonality  errors.  This  was  designated  Scenario  3  in  the  Thesis  by  A. 
Hayal  (2010).  It  was  the  basis  for  most  of  the  subsequent  analysis  of  the  impact  of  the  self¬ 
calibration  on  inertial  positioning. 


4.5.3  Results  and  Analysis 

The  effect  on  positioning  accuracy  due  to  the  IMU  self-calibration  was  analyzed  using  simulated 
GPS  outages  along  the  trajectory  of  the  cart  as  shown  in  Figure  33.  The  estimated  gyroscope  and 
accelerometer  errors  were  used  as  initial  estimates  in  a  27-state  Extended  Kalman  Filter  (EKF). 
Moreover,  standard  deviations  of  the  initial  IMU  states  were  taken  from  the  least-squares 
solution.  Positioning  errors  and  their  accuracies  were  determined  by  comparing  the  estimated 
trajectory  during  GPS  outages  to  the  true  GPS  trajectory.  Aside  from  the  original  case  which 
used  manufacturer  specified  initial  values  for  the  IMU  errors,  three  field  calibration  cases  were 
considered,  called  Scenarios  3.0,  3.3,  and  3.9.  Scenario  3.0  aimed  to  test  the  position  accuracy 
with  all  9  self  calibrated  accelerometer  errors,  only.  Therefore,  the  chosen  gyroscope  error 
parameters  were  equal  to  the  original  case  error  parameters.  Scenario  3.3  included  the  3- 
parameter  solution  for  the  gyroscope  bias  error  only,  while  the  use  of  the  9-parameter  solution 
for  all  gyroscope  errors  was  designated  as  Scenario  3.9.  The  latter  was  considered  even  though 
the  scale  factor  errors  and  non-orthogonality  errors  were  not  well  estimated  during  the  self¬ 
calibration  procedure. 

Figures  34  and  35  show  the  results  for  the  HG1700  and  the  HG1900,  respectively,  and  for 
different  lengths  of  GPS  outage.  We  found  that  it  is  possible  to  obtain  improved  position 
estimates  using  all  9  self-calibrated  accelerometer  errors.  Although  unrealistic  values  were 
obtained  for  the  gyroscope  scale  factor  errors  and  non-orthogonality  errors,  the  positioning 
accuracy  calculations  showed  that  the  full  9-parameter  gyro  error  estimates  yielded  more 
accurate  results  than  the  3-parameter  solutions  for  both  of  the  IMUs.  The  reason  for  this  has  not 
been  fully  explained,  but  it  also  occurs  only  for  the  shorter  duration  of  GPS  outage  (2s).  In 
general,  the  improvement  with  self-calibrated  errors  was  most  noticeable  for  the  HG1700,  while 
no  improvement  was  realized  (in  fact,  the  opposite  was  the  case)  for  longer  GPS  outages  being 
supported  by  the  HG1900.  This  may  be  due  to  less  stability  in  the  systematic  errors  for  the 
HG1900.  In  any  case,  for  longer  outages,  the  use  of  just  the  calibrated  gyro  biases  offers 
generally  the  best  improvement. 
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Figure  34:  HG1700  inertial  positioning  errors  with  self-calibrated  IMU  errors  according  to 
different  scenarios  (see  text)  compared  to  the  case  without  any  prior  self-calibration  prior  in  the 
field. 
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Figure  35:  HG1900  inertial  positioning  errors  with  self-calibrated  IMU  errors  according  to 
different  scenarios  (see  text)  compared  to  the  case  without  any  prior  self-calibration  prior  in  the 
field. 
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5.  Cost  analysis 


Reliable  characterization  of  MEC  requires  both  high  precision  and  high  spatial  resolution  in  the 
positioning  of  the  detector.  For  example,  Bell  (2005)  noted  that  significant  improvement  and 
efficiencies  could  be  realized  “if  a  reliable,  robust  procedure  were  available  for  precisely 
determining  the  locations  of  the  sensor  readings  without  using  the  template  [method].”  GPS 
generally  yields  1  Hz  data,  which  may  be  inadequate  for  detectors  that  pass  over  smaller  putative 
MEC  with  speeds  of  several  cm/s  or  more.  One  solution  is  to  build  a  geometrically  distributed 
array  of  detectors  that  simultaneously  passes  over  an  object  (Nelson  and  McDonald,  2001). 
Also,  laser  ranging  can  deliver  up  to  several  10s  Hz  resolution  (Leica  Geosystem1),  but  a  ranging 
system  aided  by  inertial  measurement  units  (that  provide  high  resolution,  usually  up  to  250  Hz) 
was  deemed  in  a  recent  study  (U.S.  Army  Corps  of  Engineers,  2006)  as  offering  the  most 
efficient  positioning  system  in  terms  of  resolution  and  precision.  This  synergy  of  ranging  and 
autonomous  inertial  sensor  systems  also  mitigates  various  degrading  effects  caused  by  loss  of 
line-of-sight,  multipath,  and  electromagnetic  interference. 

Our  cost  analysis  is  restricted  to  the  expense  of  determining  precise  positions  with  high 
resolution.  The  analysis  is  limited  to  tangible  startup  and  recurring  costs  and  is  based  on 
manufactured  product  costs  and  estimated  labor  costs.  There  are  numerous  more  or  less 
intangible  costs  that  can  only  be  ascertained  with  field  experience  for  any  particular  system. 
These  include  the  change  in  cost  with  a  change  in  the  rate  of  false  alarms  (which  does  not 
exclusively  depend  on  geolocation  precision),  the  differential  costs  associated  with  operating  in 
different  types  of  environment,  and  the  various  costs  of  remediation  of  different  categories  of 
MEC  in  diverse  environments. 

Another  important  intangible  aspect  of  the  analysis  is  the  operational  complexity  associated 
with  a  particular  type  of  terrain  to  be  surveyed.  All  geolocation  systems  depend  on  an  external 
reference;  that  is,  no  system  (of  reasonable  technology  and  cost)  can  operate  autonomously  with 
the  kind  of  accuracy  needed  for  initial  screening,  mapping,  and  discrimination  and 
characterization  of  MEC.  The  type  of  terrain  and  environment  dictates  the  accessibility  and, 
more  importantly,  the  degree  of  sustained  access  to  such  external  references.  As  in  many  cases 
of  comparative  analyses  of  different  system  there  are  tradeoffs  to  which  costs  may  or  may  not  be 
particularly  sensitive.  Ultimately,  we  argue  that  for  any  geolocation  system  the  inclusion  of  an 
IMU  is  a  cost-effective  enhancement,  if  only  a  safeguard,  for  precision  positioning  of  MEC 
detectors. 

We  compare  three  competing  systems  that  offer  high-resolution  positioning  in  a  variety  of 
possibly  GPS-degraded  environments.  The  system  that  is  the  basis  of  our  analysis  for  this 
project  is  the  integrated  GPS/IMU  geolocation  system,  comprising  a  geodetic  quality  GPS 
receiver  (and  antenna)  and  a  medium-grade  IMU  (such  as  the  Honeywell  HG1700  or  HG1900). 
Another  external  positioning  system  is  the  laser-based  system  whereby  the  roving  (detection) 
system  is  positioned  using  trilateration  and  triangulation  with  respect  to  a  number  of  fixed 
ground  stations  occupied  by  laser-distance  measuring  devices  (robotic  total  stations).  The  third 
system  is  an  acoustic  system,  where  distances  are  measured  between  an  ultrasonic  signal 
transmitter  on  the  geophysical  sensor  unit  and  several  stationary  receivers  in  the  field.  In 
essence,  all  external  positioning  systems  can  be  thought  of  as  a  number  of  benchmarks  with 
known  positions  transmitting  (receiving)  a  signal  to  (from)  the  roving  target,  which  either 


1  http://www.leica-geosystems.com/ 
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through  active  or  passive  participation  leads  to  a  determination  of  its  position  in  a  well  defined 
coordinate  system.  With  the  ground  laser  trilateration  system,  the  benchmarks  are  stations  in  the 
field  and  the  roving  detection  system  passively  reflects  the  transmitted  laser  pulses;  with  GPS, 
the  benchmarks  are  the  satellites  that  transmit  their  position  and  time  continuously  and  the  roving 
system  receives  and  utilizes  these  signals  to  determine  its  position. 

In  all  cases,  the  line  of  sight  between  the  benchmark  (satellite)  and  the  roving  detection 
system  must  be  unobstructed  for  the  particular  wavelengths  of  the  electro-magnetic  radiation  on 
which  the  positioning  system  is  based.  Therefore,  in  principle,  the  IMU  could  serve  as  an  aid  to 
each  ranging  system  in  the  event  that  the  line-of-sight  (or  acoustic  range)  is  temporarily 
obstructed,  interrupted,  or  degraded  because  of  various  environmental  factors.  However,  for 
comparison  purposes,  we  assume  that  the  GPS  degradation  is  the  most  common  and  serious, 
since  the  other  ranging  systems  are  ad  hoc  systems  designed  specifically  for  the  application  in  a 
particular  environment.  GPS,  on  the  hand,  would  be  considered  a  viable  system  in  any 
environment  that  yields  at  least  a  reasonable  access  to  the  satellite  signals.  In  that  sense  an 
integrated  GPS/IMU  may  be  compared  to  the  other  stand-alone  ranging/geolocation  systems. 

Table  7  summarizes  the  basic  features  of  three  representative  geolocation  systems  as 
proposed  or  realized  as  described  in  the  final  report  by  the  US  Army  Corps  of  Engineers  (2006), 
and  as  studied  by  us.  All  three  systems  rely  on  GPS  or  equivalent  to  perform  absolute 
positioning.  The  principal  innovation  is  in  precision  relative  position  of  the  UXO  detection 
system.  The  systems  differ  in  operational  principle,  where  the  roving  unit  for  the  GPS/IMU  and 
the  ultrasonic  systems  contain  the  main  geolocation  sensors,  while  for  the  robotic  laser-tracking 
system  the  principal  equipment  is  a  stationary  tracking  device.  As  such  the  costs  are  not 
comparable  in  a  direct  way.  While  the  GPS/IMU  system  is  completely  mobile,  the  laser  tracking 
system  is  restricted  to  the  range  of  the  total  station,  until  it  can  or  must  be  moved  (which  requires 
additional  manpower).  Similarly,  other  local  ranging  systems,  exemplified  by  the  ultrasonic 
system  require  repositioning  of  equipment  at  reference  station  as  defined  by  the  survey  area. 


Table  7:  Cost  com] 

parison  of  geolocation  systems 

System 

Component  equipment 

Cost  (ROM 
estimates) 

notes 

GPS/IMU 

GPS  receiver  and  antenna 
Inertial  Measurement  Unit 
(HG1900) 

$29000 

$10000 

ancillary  equipment  estimated 
cost  of  about  $500 

Robotic  Total 
Station  (RTS) 

Laser  total  station,  plus 
robotic  tracking, 
LeicaTSPllOO 

$35000 

cost  for  single  tracking  station; 
excludes  cost  of  absolute 
positioning  equipment 

GPS/ultrasonic 

ranging 

GPS  receiver  and  antenna 
Ultrasonic  units 

$29000 

$1500 

ancillary  equipment  estimated 
cost  of  about  $6000 

The  principal  disadvantage  of  the  IMU  system  for  longer-term  geolocation  is  the  inherent 
drift  in  the  navigation  solution,  which  is  absent  in  all  ranging  systems.  For  short-term 
interpolation  between  reference  points  (say,  less  10  s),  the  medium-grade  IMU  with  appropriate 
post-processing  performs  at  levels  about  an  order  of  magnitude  worse  than  the  best  ranging 
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(laser)  systems.  The  main  advantage  of  the  IMU  is  its  complete  autonomy;  requiring  no 
additional  base  stations  and  being  unaffected  by  most,  if  not  all  environmental  conditions. 

We  showed  in  our  various  laboratory  and  field  tests  that  the  HG1700  is  the  slightly  more 
precise  than  the  HG1900.  However,  the  cost  is  60%  greater  and  may  not  be  warranted  under 
most  scenarios.  Considering  the  cost  of  the  medium-grade  IMU,  HG1900,  relative  to  the  typical 
cost  of  an  alternative  geolocation  system,  one  may  argue  that  even  with  the  alternative  systems, 
the  rover  should  include  an  IMU  to  bridge  any  possible  outages  in  the  ranging  system  (although 
most  systems  claim  almost  full-time  integrity,  any  ranging  system  is  always  susceptible  to 
obstruction  or  interference).  The  data  processing  techniques  developed  under  this  project  are 
applicable  to  the  integration  of  IMU  data  with  any  ranging  system. 
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6.  Summary  and  Conclusions 


To  satisfy  precise  positioning  requirements  for  MEC  detection  and  characterization,  we 
considered  the  integrated  ranging/IMU  geolocation  system.  Primarily,  we  studied  the  algorithms 
that  optimally  combine  IMU  data  with  updates  provided  by  GPS.  Our  analysis  proceeded  from 
simulations  to  field  data  and  included  various  novel  estimation  tools  developed  recently  for 
inertial  navigation  of  land  and  airborne  vehicles.  We  also  studied  the  static  self-calibration  of 
systematic  IMU  errors  in  the  field.  The  new  estimation  techniques  focused  on  the  non-linearity 
of  the  trajectory  dynamics  and  the  possible  non-Gaussianity  of  the  random  instrument  errors. 
Accommodating  these  characteristics  requires  more  general  Bayesian  estimation  than  developed 
for  the  extended  Kalman  filter.  Based  on  simulated  trajectories  typical  of  MEC  ground  surveys, 
as  well  as  cart-based  and  hand-held  trajectories  in  the  laboratory  and  actual  MEC  detection 
trajectories  performed  in  the  field,  we  analyzed  the  unscented  Kalman  filter,  the  unscented 
particle  filter,  a  hybridization  of  these  non-linear  filters  and  the  extended  Kalman  filter,  known  as 
the  Rao-Blackwellized  filter,  and  various  modifications,  including  adaptive  error  techniques, 
neural  network  applications,  and  wave-correlation  filters  (for  the  case  the  dual  IMUs  are 
utilized).  We  also  considered  the  smoothing  version  of  each  of  these  filters,  wherein  the 
estimated  IMU  trajectory  is  controlled  optimally  over  the  entire  inter-update  interval  by  the 
updates  at  its  endpoints. 

Our  simulations  showed  that  the  unscented  Kalman  filter  (based  on  the  unscented 
transformation  that  bypasses  a  linearization  of  the  error  state  dynamics  and  of  the  observation 
updates)  performs  consistently  better  than  the  standard  extended  Kalman  filter,  particularly  along 
curved  trajectories.  These  improvements  in  filter  strategy  were  realized  especially  when  the 
interval  of  the  ranging  solution  update  was  several  seconds  (simulating  an  outage  due  to  signal 
occlusions)  and  if  the  ranging  solution  was  degraded  (simulating  various  possible  causes). 

Using  particle  filters  avoids  the  Gaussianity  assumption  and  our  tests  showed  that  these 
filters  are  particularly  useful  if  the  driving  noise  of  the  system  has  an  asymmetric  distribution.  In 
this  case,  the  UKF  and  EKF  perform  comparably,  but  the  UPF  yields  significantly  improved 
position  accuracy.  The  UPF  results  were  generally  insensitive  to  the  number  of  particles,  and 
improvement  could  be  obtained  (in  the  case  of  longer  GPS  outages)  using  adaptive  techniques 
that  compensated  for  the  UT’s  assumption  of  symmetric  noise  probability  densities. 

From  the  various  filters  tested,  we  find  that  achieving  few  centimeters  of  positioning 
accuracy  in  dynamic  environments  requires  non-linear  filters,  such  as  the  UKF  or  UPF.  These 
filters  cannot  overcome  the  natural  accumulation  of  IMU  errors  as  the  ranging  solution  update 
interval  increases.  However,  in  all  cases  the  new  non-linear  filters  performed  better  than  the 
standard  EKF.  The  best  performance  among  all  filters  tested  was  obtained  by  the  AUPF  which 
accommodates  non- symmetric  sensor  errors  as  well  as  highly  dynamic  trajectories. 

These  results  were  validated  using  actual  cart-based  and  hand-held  IMU  trajectories  in  the 
laboratory  that  simulate  typical  field  trajectories.  Laboratory  tests  were  conducted  using  IMUs 
with  different  capabilities  and  along  trajectories  with  different  dynamics.  The  positioning 
performance  was  evaluated  for  linear  and  non-linear,  as  well  as  adaptive  and  neural- network- 
aided  adaptive  Kalman  filter/smoothers.  As  expected,  the  nonlinear  smoothers,  developed  from 
a  combination  of  adaptive  unscented  Kalman  filter  and  RTS  smoother  (AUKS)  yielded  superior 
performance  over  the  standard  adaptive  extended  Kalman  smoother.  The  neural- network  aiding, 
in  particular,  tended  to  decrease  the  difference  in  performance  between  benign  and  dynamic 
components  of  the  trajectory. 
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A  Handheld  Geolocation  System  (HGS)  using  a  tactical-grade  IMU  was  designed  to  obtain 
precise  positions  of  a  UXO  detection  system  applied  in  relatively  small  areas  (1.0m  by  1.0m). 
Since  the  test  was  operated  in  a  closed  environment  where  no  GPS  signal  was  available,  an 
automatic  position  tracking  system  was  designed  and  implemented.  This  would  be  comparable 
to  GPS  if  a  GPS  antenna  were  mounted  on  top  of  the  handheld  system.  As  before,  the  UPF 
performed  generally  best  among  all  other  (linear  and  nonlinear)  filters.  However  since  the  UPF 
needs  a  large  number  of  samples  to  represent  the  a  posteriori  statistics  in  high-dimensional  space, 
the  Rao-Blackwellized  Particle  Filter  (RBPF)  was  tested  as  an  alternative  to  increase  the 
efficiency  of  the  particle  filter.  The  corresponding  nonlinear  smoothers  that  are  based  on 
forward/backward  filtering  techniques  (EKS,  UPS,  and  RBPS)  were  also  tested  and  analyzed  for 
two  typical  local  handheld  detection  platform  trajectories  (sweep  and  swing).  On  the  whole, 
position  accuracy  improvements  were  achieved  by  applying  nonlinear  filter-based  smoothing 
techniques  (UPS  and  RBPS)  in  both  the  straight  and  curved  sections  of  the  trajectories.  The 
handheld  geolocation  system  with  a  nonlinear  filter-based  smoother  achieved  the 
characterization  and  discrimination  level  of  accuracy  if  the  update  rate  of  control  points  is  less 
than  0.5Hz  and  1Hz  for  the  sweep  and  swing  modes  respectively.  Although  the  data  collection 
using  the  swing  operation  can  be  done  in  shorter  time,  the  sweep  operation  is  generally  better 
than  the  swing  because  the  dynamics  of  swing  operation  is  highly  dynamic. 

The  improved  performance  of  the  IMU/GPS  geolocation  system  using  the  novel  estimation 
methods  was  again  verified  in  the  field  on  towed-vehicle  and  cart-based  MEC  detection 
platforms.  In  collaboration  with  NRL  personnel,  we  conducted  tests  at  NRL’s  UXO  test  site 
located  at  the  Army  Research  Laboratory  Blossom  Point  Facility  in  Maryland.  In  addition,  to  the 
non-linear  filters  and  smoothers,  we  employed  the  wave -correlation  filter  (which,  however, 
requires  a  dual  IMU  system),  as  well  as  simple  endpoint  matching  as  an  alternative  to  the 
smoothing  algorithms.  The  wave-correlation  filter  provided  only  modest  improvement,  while 
the  endpoint  matching  yielded  similar  performance  to  the  smoothing  algorithms,  due  to  the 
relatively  benign  (straight-line)  trajectories  that  were  tested. 

Finally,  we  tested  a  self- calibration  method  in  the  field  using  our  own  cart-based  geolocation 
system  that  included  a  dual  IMU  and  GPS  receiver  and  antenna.  The  systematic  errors  of  the 
inertial  sensors  could  be  calibrated  in  the  static  mode  by  orienting  their  sensitive  axes  in  various 
directions  relative  to  the  Earth’s  gravity  and  rotation  vectors.  In  this  way  biases,  scale  factor 
errors,  and  non-orthogonality  errors  in  the  accelerometer  outputs  were  estimated,  as  well  as  bias 
errors  in  the  gyro  outputs.  It  was  also  demonstrated  that  these  prior  calibrations  in  the  field 
improved  the  inertial  trajectory  performance,  at  least  for  short  (few  second)  update  intervals. 

A  rudimentary  cost  analysis  of  the  IMU/GPS  geolocation  system  showed  that  it  is 
competitive  in  terms  of  cost  with  respect  to  alternative  systems  that  attempt  to  improve  or  aid  the 
standard  GPS  ranging  system.  However,  the  primary  function  of  IMU  aiding  is  more  tuned  to 
short-term  bridging  of  GPS  outages  and  providing  much  higher  resolution  than  GPS  can 
typically  offer,  similar  to  or  better  than  laser-tracking  systems.  IMU  aiding  is  also  distinguished 
by  the  fact  that  it  requires  no  external  support  system  and  thus  is  completely  autonomous  and 
totally  immune  to  external  environmental  influences. 
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